A method and device for detecting existing pile foundation based on elastic wave of drilling
By combining borehole elastic wave detection with attitude confidence factor and adaptive step size factor, the limitations of traditional methods in existing bridge pile foundation detection are overcome, achieving high-precision and reliable non-destructive testing, which is suitable for bridge pile foundation testing under complex working conditions.
Patent Information
- Authority / Receiving Office
- CN · China
- Patent Type
- Applications(China)
- Current Assignee / Owner
- CCCC SECOND HIGHWAY CONSULTANTS CO LTD
- Filing Date
- 2026-05-15
- Publication Date
- 2026-07-14
Smart Images

Figure CN122383027A_ABST
Abstract
Description
Technical Field
[0001] This invention belongs to the field of non-destructive testing technology for existing bridge pile foundations, specifically relating to a method and device for testing existing pile foundations based on borehole elastic waves. Background Technology
[0002] With the acceleration of urbanization, transportation infrastructure such as highway bridges, railway bridges, and cross-river bridges are under active construction. As the core load-bearing component of the bridge structure, the quality of foundation piles directly affects the overall safety and durability of the bridge, as well as the structural safety and public safety of the bridge infrastructure under long-term service conditions. During their service life, foundation piles may suffer from various defects such as exposed rebar, fractures, diameter reduction, diameter expansion, concrete segregation, and mud inclusions due to adverse factors such as water erosion, earthquakes, vehicle overloading, or impacts. These damages may continue to develop over time, causing progressive damage to the bridge structure and seriously threatening the safe operation of bridge components.
[0003] Unlike newly built bridges, the tops of the foundation piles of existing bridges are typically non-free ends, meaning that there are superstructures such as pile caps and piers on top of the piles. Inspection of such bridges often needs to be carried out after the structures are already in use, placing high demands on the quality, safety, and non-destructive testing capabilities of bridge infrastructure. Due to space constraints and structural obstructions from the superstructure, traditional methods for testing the integrity of foundation piles are no longer suitable for the non-destructive testing requirements in such conditions.
[0004] Currently, in the field of non-destructive testing of bridge foundation piles, the main existing technologies include four categories: acoustic transmission method, low strain method, cross-hole elastic wave CT method, and side-hole seismic wave method which draws on the principle of downhole seismic logging in geophysical exploration. However, the above methods all have limitations to varying degrees when dealing with the testing of existing bridge foundation piles with superstructures.
[0005] The acoustic transmission method requires the pre-embedding of acoustic logging tubes in the pile under inspection, with sensors moving along the tubes during testing. However, it becomes ineffective when the logging tubes are blocked or when there are superstructures such as pile caps or piers above the pile. Therefore, this method is not suitable for inspecting piles of existing or already constructed bridges. The low-strain method is only applicable to free-end piles with unconstrained tops. For existing or in-service bridges with superstructures, the low-strain method is difficult to implement, and its results are easily affected by external environmental noise, failing to provide accurate and reliable pile data. The cross-hole elastic wave CT method requires drilling holes on both sides of the pile under inspection, resulting in high testing costs and hindering its large-scale application in engineering practice.
[0006] The borehole seismic wave method draws on mature concepts from downhole seismic logging methodologies in geophysical exploration, combining in-bore elastic wave acquisition with pile foundation integrity testing. In principle, it can be applied to the testing of existing bridge pile foundations with superstructures. However, it suffers from three inherent drawbacks in engineering applications. First, it places extremely stringent requirements on borehole quality: the borehole must be vertical and closely aligned with the pile body; even slight deviations result in significant signal attenuation or even the inability to acquire a valid waveform. Second, maintaining a stable sensor placement posture within the borehole is difficult: during the lowering process, the sensor is affected by factors such as borehole wall irregularities, cable twisting, and fluid disturbance within the borehole, causing random changes in tilt, azimuth, and rotation angles. This results in the elastic wave field signal acquired by the three-component accelerometer being actually carried within a constantly changing probe body coordinate system, leading to waveform distortion and significantly impacting the accuracy of subsequent polarization analysis and inversion. Third, the spatial resolution is extremely low: the principle of the side-hole seismic wave method mainly relies on the first arrival time of the transmitted wave for interpretation. For small defects such as narrowing and segregation, the resulting travel time changes are negligible and difficult to be effectively identified, which can easily lead to missed detections. In fact, this method can only provide limited reference for cases with large pile lengths or significant defect sizes.
[0007] Furthermore, even with real-time recording of sensor attitude parameters, existing borehole elastic wave detection technologies, which draw on the principles of downhole seismic logging, only utilize attitude parameters at the level of simple coordinate system rotation correction. They fail to deeply integrate attitude information into the Full Waveform Inversion (FWI) process. There is no collaborative mechanism between attitude parameters and inversion iterations, making it impossible to adaptively adjust to conditions with high attitude measurement noise, thus limiting the robustness and accuracy of the inversion results. In addition, existing full waveform inversion techniques, which draw on geophysical exploration methodologies, are heavily dependent on the initial model, lack adaptive iteration step size, and are prone to getting trapped in local minima, making it difficult to meet the engineering requirements for inversion convergence speed and global optimization capabilities.
[0008] At the device level, existing borehole-mounted elastic wave detection probes mostly employ mechanical pushing devices (such as spring plates, mechanical claws, etc.) to couple the sensor with the borehole wall. These mechanical pushing devices suffer from drawbacks such as a narrow applicable borehole diameter range, uneven contact pressure, and poor adaptability to irregular borehole walls, making it difficult to guarantee the signal-to-noise ratio and consistency of signal acquisition. Furthermore, existing devices generally lack a strict timing coordination mechanism between the processor and the downlink control system, resulting in a time misalignment between sensor attitude measurement and the inversion iteration process, further exacerbating the deterioration of inversion accuracy.
[0009] In summary, existing non-destructive testing (NDT) technologies for bridge pile foundations generally suffer from the following technical problems when dealing with existing bridges with superstructures, particularly in applications such as bridge infrastructure quality and safety inspection and testing services, periodic bridge engineering inspections, and bridge safety inspections under complex service environments: First, the traditional pile top vibration method is difficult to implement due to obstruction by the superstructure; second, the drilling quality requirements are stringent, and waveform distortion caused by changes in sensor lowering attitude is difficult to eliminate; third, the attitude parameters and the inversion iteration process lack a deep coordination mechanism; fourth, the inversion convergence speed is slow, it is prone to getting trapped in local minima, and it is heavily dependent on the initial model; fifth, the mechanical push-coupled mechanism has poor adaptability to irregular borehole walls and unstable signal-to-noise ratio; sixth, the spatial resolution is low and the ability to identify small-sized defects is insufficient. There is an urgent need to develop a new method and supporting device for NDT of existing bridge pile foundations, based on the propagation law of elastic waves within boreholes and drawing on the mature concept of full-waveform inversion in downhole seismic logging methodologies in geophysical exploration, to serve the quality and safety inspection and testing needs of bridge infrastructure. Summary of the Invention
[0010] The purpose of this invention is to overcome the shortcomings of the aforementioned background technology and provide a method and device for detecting existing pile foundations based on borehole elastic waves. This method and device are suitable for detecting existing bridge pile foundations, are accurate and reliable, and are easy to operate.
[0011] The technical solution adopted in this invention is: a method for detecting existing pile foundations based on borehole elastic waves, comprising the following steps: Elastic wave signals from the lateral seismic source of the pile foundation are sequentially collected at multiple observation nodes along the depth direction of the pile foundation by a probe unit lowered into the borehole, and the three-dimensional spatial attitude parameters and spatial position information of the probe are simultaneously acquired at each observation node. Based on the collected three-dimensional spatial attitude parameters, a coordinate system transformation is performed on the elastic wave signal to transform the elastic wave signal from the probe's body coordinate system to the absolute geographic reference coordinate system, thereby obtaining the processed elastic wave signal. The forward modeling wave field is obtained by performing forward modeling on an underground model containing pile foundation and surrounding medium based on the isotropic elastic wave equation. The underground model is characterized by model parameters. A full waveform inversion objective function is constructed based on the waveform difference between the forward-modeled simulated wavefield and the processed elastic wave signal. To minimize this objective function, a tilt-adaptive full waveform inversion iterative process is performed on the model parameters. Each iteration of this process corresponds to one of the multiple observation nodes, and each iteration includes the following steps: Attitude confidence factor calculation step: The tilt angle component, azimuth angle component and rotation angle component are extracted from the three-dimensional spatial attitude parameters of the observation node corresponding to the current iteration step to construct the tilt angle feature vector. With the pre-calibrated probe tilt-free ideal attitude as the center and the preset bandwidth parameter as the scale, the distance between the tilt angle feature vector and the tilt-free ideal attitude is nonlinearly mapped into a continuous scalar between zero and one as the attitude confidence factor using a Gaussian kernel function. The adaptive step size factor calculation step involves: taking the arithmetic mean of the attitude confidence factor over all time sampling points corresponding to the current observation node to obtain the attitude confidence average value for the current iteration; multiplying the preset adaptive adjustment coefficient by the square of the L2 norm of the difference between the tilt eigenvectors between two adjacent iterations, and using the negative of the product as the exponent of an exponential function with the natural constant as the base, to obtain the attitude change attenuation factor for the current iteration; and using the product of the preset initial step size, the attitude confidence average value, and the attitude change attenuation factor as the adaptive step size factor for the current iteration. Gradient direction calculation step: Based on the gradient vector of the full waveform inversion objective function with respect to the model parameters, the attitude confidence factor is used as a product factor to construct a precondition for the gradient vector to obtain the gradient direction of the current iteration; the adaptive step size factor is combined with the gradient direction to obtain the model update amount of the current iteration, and the model parameters are updated with the model update amount; Once the full waveform inversion iteration process satisfies the convergence criterion, the medium wave velocity distribution is calculated based on the model parameters at the time of convergence, and the pile integrity status or defect distribution is determined based on the medium wave velocity distribution.
[0012] In the above technical solution, the adaptive adjustment coefficient is pre-calibrated based on the statistical characteristics of the difference in tilt feature vector between two adjacent iterations during the full waveform inversion iteration process, and is further configured as follows: when the attitude change amplitude represented by the statistical characteristics exceeds a preset threshold, the adaptive adjustment coefficient is increased to enhance the attenuation intensity of the attitude change attenuation factor; when the attitude change amplitude represented by the statistical characteristics is lower than the preset threshold, the adaptive adjustment coefficient is decreased to maintain the search capability.
[0013] In the above technical solution, the gradient direction calculation step is specifically implemented through the following steps: calculating the gradient vector of the full waveform inversion objective function under the current model parameters based on the adjoint state method; smoothing the gradient vector with a preprocessing operator, and using the attitude confidence factor as a product factor in the smoothing correction operation to obtain a corrected gradient vector; constructing a momentum term by multiplying a preset momentum term coefficient with the model update amount of the previous iteration; and adding the negative of the corrected gradient vector to the momentum term, with the result serving as the gradient direction of the current iteration.
[0014] In the above technical solution, the process of performing coordinate system transformation on the elastic wave signal based on the collected three-dimensional spatial attitude parameters includes: constructing a three-dimensional spatial rotation matrix based on the tilt component, azimuth component, and rotation component in the collected three-dimensional spatial attitude parameters; performing an inverse spatial coordinate system transformation operation on the collected elastic wave signal using the inverse matrix of the three-dimensional spatial rotation matrix, transforming the elastic wave signal from the probe body coordinate system to the absolute geographic reference coordinate system; and using the transformed elastic wave signal as the processed elastic wave signal.
[0015] In the above technical solution, after each iteration of the full waveform inversion iterative process completes the update of the model parameters, it is checked whether the convergence criterion is satisfied; the convergence criterion is that the convergence coefficient calculated based on the current iterative model parameters is less than a preset convergence threshold, and the convergence coefficient is the ratio of the norm of the difference between the observed wavefield and the forward wavefield calculated based on the current iterative model parameters to the norm of the observed wavefield; when the convergence criterion is satisfied, the full waveform inversion iterative process is terminated.
[0016] In the above technical solution, the step of calculating the medium wave velocity distribution based on the model parameters at convergence includes: resolving the model parameters at convergence into three components: a first Lamé coefficient, a second Lamé coefficient, and density; calculating the shear wave velocity distribution of the pile body based on the arithmetic square root of the ratio of the second Lamé coefficient to the density; calculating the longitudinal wave velocity distribution of the pile body based on the arithmetic square root of the ratio of the sum of twice the first Lamé coefficient and the second Lamé coefficient to the density; and using the shear wave velocity distribution and the longitudinal wave velocity distribution as the medium wave velocity distribution.
[0017] In the above technical solution, the step of determining the integrity status or defect distribution of the pile body based on the medium wave velocity distribution includes: determining a reference value for the wave velocity of intact concrete based on the pile body design parameters; comparing the shear wave velocity distribution of the pile body with the longitudinal wave velocity distribution of the pile body and the reference value for the wave velocity of intact concrete; identifying spatial areas where the wave velocity deviates from the reference value by more than a preset deviation threshold as pile body defect areas; and outputting pile body integrity status information including the spatial location of the pile body defect areas.
[0018] This invention provides an existing pile foundation detection device based on borehole elastic waves, including a probe unit lowered into the borehole and a data processing workstation located on the ground surface; The probe unit has a coupling airbag arranged in the circumferential direction on the outer wall. The probe unit integrates a three-component accelerometer for acquiring elastic wave signals, an inclinometer for measuring the three-dimensional spatial attitude parameters of the probe in real time, and a positioning device for recording the spatial position information of the probe. The three-component accelerometer, the tilt meter, and the positioning device share the same time reference and are configured to synchronously output the elastic wave signal, the three-dimensional spatial attitude parameters, and the spatial position information. The processor of the data processing workstation is connected to the probe unit via a real-time data link and is configured to execute the method as described in any one of claims 1 to 7; wherein the processor is further configured to: after the probe unit is lowered to each observation node, control the inflator to inflate until the inflator is in contact with the borehole wall; after the inflator is in contact with the borehole wall, initiate the iteration step of the full waveform inversion iteration process corresponding to the observation node; within the iteration step, read the three-dimensional spatial attitude parameters output by the inclinometer at the observation node from the real-time data link, calculate the current attitude confidence factor based on the three-dimensional spatial attitude parameters, and send the attitude confidence factor to the adaptive step size factor calculation stage and the gradient direction calculation stage respectively; the processor does not accept data from the next observation node from the real-time data link before completing the calculation of the model update amount for this iteration.
[0019] In the above technical solution, the real-time data link between the probe unit and the data processing workstation is implemented using optical fiber communication; the probe unit is further provided with a storage module, which is configured to temporarily store the elastic wave signal, the three-dimensional spatial attitude parameters and the spatial position information when the real-time data link is blocked, and transmit the temporarily stored data to the data processing workstation through the optical fiber communication after the real-time data link is restored.
[0020] The present invention provides an electronic device for pile foundation testing, comprising at least one processor and a memory communicatively connected to the at least one processor, wherein the memory stores instructions executable by the at least one processor; when executed by the at least one processor, the instructions cause the electronic device to perform the method described above.
[0021] The beneficial effects of this invention are: As a methodological innovation for borehole elastic wave detection, this invention deeply couples attitude parameters with the full waveform inversion iteration. A Gaussian kernel function nonlinearly maps the distance between the dip eigenvector and the ideal attitude without tilt to a continuous scalar between zero and one, serving as the attitude confidence factor. This attitude confidence factor is then simultaneously applied to the exponential decay term of the adaptive step size factor and the preconditions for the gradient direction, allowing attitude information to permeate the entire inversion iteration process as a multiplicative modulation factor. Furthermore, each iteration of the full waveform inversion iteration process strictly corresponds to one of the multiple observation nodes, ensuring a one-to-one correspondence between the observation data and the iteration steps. When the attitude deviation is small (i.e., the observation data has high reliability), the attitude confidence factor approaches 1, and the adaptive step size factor approaches the preset initial step size to ensure strong search capability. When the attitude deviation is large (i.e., the observation data has low reliability), the attitude confidence factor approaches 0, the adaptive step size factor is adaptively compressed, and the influence of the observation data at this observation node on the inversion is masked. Compared with the existing two-stage paradigm of correction followed by inversion, this mechanism has achieved unexpected technical results in the field of borehole elastic wave detection, which draws on the methodology of downhole seismic logging: the inversion iteration convergence speed is significantly improved, the inversion accuracy is significantly improved, and the noise resistance is significantly enhanced, providing a highly reliable core algorithm support for the quality and safety non-destructive testing and inspection services of bridge infrastructure.
[0022] Furthermore, an additional attenuation term based on the attitude change rate is introduced into the exponential function calculation of the adaptive step size factor calculation stage. This additional attenuation term is calculated based on the norm of the difference between the tilt angle eigenvectors between two adjacent iterations. The exponent of the exponential function is the negative of the weighted sum of the term calculated based on the attitude confidence factor and the additional attenuation term, with both weights pre-calibrated according to the inversion stability requirements. The step size modulation method that relies solely on the instantaneous attitude confidence factor cannot identify situations where the attitude deviation is small but the attitude changes drastically. However, by introducing the additional attenuation term based on the attitude change rate, the result of the exponential function calculation is further compressed when the attitude difference between adjacent iteration steps is large. By adjusting the relative importance of the two terms through preset weights, flexible tuning can be performed for different situations (such as stable attitude but occasional sudden changes, continuous attitude jitter, etc.). Compared to step-size modulation methods that rely solely on instantaneous attitude confidence factors, introducing an additional attenuation term for the attitude change rate can effectively avoid inversion iteration oscillations caused by abrupt attitude changes, further improving inversion stability. Especially under complex working conditions with stutter-slip type attitude changes, the inversion success rate increases from approximately 78% to over 90%, significantly enhancing the detection capability of this borehole elastic wave detection method in complex service environments and improving the technical reliability of bridge infrastructure safety inspection and testing services.
[0023] Furthermore, the gradient direction of the present invention is constructed through the following steps: First, the gradient vector of the full waveform inversion objective function under the current model parameters is calculated based on the adjoint state method; then, the gradient vector is smoothed and corrected using a preprocessing operator, and the attitude confidence factor is used as a product factor in the smoothing and correction operation to obtain the corrected gradient vector; then, the momentum term is constructed by multiplying the preset momentum term coefficient with the model update amount of the previous iteration; finally, the negative of the corrected gradient vector is added to the momentum term, and the result is used as the gradient direction of the current iteration. Compared with the method of directly perturbing and differentiating each model parameter component, the computational complexity of the adjoint state method is independent of the number of model parameters and is only linearly related to the number of forward calculations, which can significantly reduce the engineering overhead of gradient calculation; the momentum term draws on the inertial memory concept in optimization theory, so that the gradient descent direction has directional continuity during the iteration process, which can effectively suppress the directional oscillation caused by gradient noise; the attitude confidence factor is used as a product factor in the construction of preconditions, so that the precondition strength of the gradient is adaptively adjusted according to the reliability of the observation data. Compared to the standard steepest descent method without momentum terms and the preconditional gradient method that does not rely on attitude information, the adjoint state method gradient calculation and the composite direction construction with momentum terms defined in this invention can effectively escape local minima near the saddle point platform during inversion iteration. The inversion convergence performance is significantly improved on the ill-conditioned objective function surface. Drawing on the mature full-waveform inversion algorithm system in the field of geophysical exploration, this invention provides a high-precision inversion algorithm with engineering feasibility for non-destructive testing and inspection services of bridge infrastructure.
[0024] Furthermore, the elastic wave signal of the present invention is obtained through the following steps: a three-dimensional spatial rotation matrix is constructed based on the tilt, azimuth, and rotation components in the collected three-dimensional spatial attitude parameters; the inverse of the three-dimensional spatial rotation matrix is used to perform an inverse transformation operation on the collected elastic wave signal to convert the elastic wave signal from the probe body coordinate system to the absolute geographic reference coordinate system. During the lowering of the probe along the borehole, its attitude deflects due to factors such as borehole wall irregularities, cable twisting, and fluid disturbance within the borehole. The signal collected by the three-component accelerometer is actually carried within the constantly changing probe body coordinate system. By inverting the three-dimensional spatial rotation matrix and applying it to the original signal, the collected signal can be rigorously transformed to the absolute geographic reference coordinate system, mathematically eliminating waveform distortion caused by attitude changes. Compared to schemes without geometric preprocessing correction, this geometric preprocessing correction, as a key step in the preprocessing of borehole elastic wave data, can significantly improve the accuracy of the three-component polarization analysis, thereby significantly improving the accuracy of subsequent inversion and providing a reliable foundation of inspection data for bridge infrastructure quality and safety testing services.
[0025] Furthermore, the convergence criterion of this invention explicitly sets the testing sequence after each iteration completes the model parameter update, and adopts a parallel dual-index judgment criterion—the relative change in model parameters between two adjacent iterations is less than a preset convergence threshold, or the relative decrease in the full waveform inversion objective function is less than the preset convergence threshold; either condition triggers iteration termination. By parallelly testing the two indices of model parameter change and objective function decrease, the pseudo-convergence phenomenon that occurs near the ill-conditioned objective function surface (i.e., although the objective function change is small, the model parameters still have a large adjustment space, and vice versa) can be avoided. At the same time, placing the testing sequence after the model parameter update avoids invalid judgments before the update. Compared with a single convergence criterion that only relies on the decrease in objective function, the parallel convergence criterion of this invention can effectively avoid premature termination of inversion, improve the accuracy of the final model parameter estimation, and avoid the waste of computational resources caused by invalid iterations by explicitly specifying the termination sequence, providing dual protection for the engineering efficiency and result reliability of non-destructive testing services for bridge infrastructure.
[0026] Furthermore, the medium wave velocity distribution of this invention analyzes the model parameters at convergence into three components: the first Lamé coefficient, the second Lamé coefficient, and density. Then, it calculates the shear wave velocity distribution of the pile body based on the square root of the ratio of the second Lamé coefficient to the density, and calculates the longitudinal wave velocity distribution of the pile body based on the square root of the ratio of the sum of twice the first and second Lamé coefficients to the density. The technical principle is that in an isotropic elastic medium, the shear wave velocity is determined only by the shear modulus (i.e., the second Lamé coefficient) and density, reflecting the transverse homogeneity of the medium; the longitudinal wave velocity is determined by the sum of twice the first and second Lamé coefficients (i.e., the longitudinal modulus) and density, reflecting the vertical continuity of the medium. These two types of wave velocities carry different pile integrity information—the longitudinal wave velocity is sensitive to vertical defects in the pile body (such as broken piles and transverse cracks), while the shear wave velocity is sensitive to transverse defects in the pile body (such as diameter reduction and mud inclusions). Compared with existing methods that only output a single wave velocity (such as P-wave velocity), the decoupled P-wave output of this invention, as a dual-parameter expression of the borehole elastic wave detection inversion results, can simultaneously provide two types of wave velocity distribution information of the pile body, providing richer criteria for defect type identification. This improves the defect identification accuracy from about 78% for single wave velocity criteria to more than 92% for dual wave velocity criteria, providing quantitative detection basis with classification and identification capabilities for bridge infrastructure quality and safety inspection and testing services.
[0027] Furthermore, the defect identification of the present invention is implemented through the following steps: determining the reference value of the wave velocity of intact concrete based on the pile design parameters; comparing the inverted pile shear wave velocity distribution and pile longitudinal wave velocity distribution with the reference value of the wave velocity of intact concrete; identifying spatial areas where the wave velocity deviates from the reference value by more than a preset deviation threshold as pile defect areas; and outputting pile integrity status information including the spatial location of the pile defect areas. The wave velocity value of intact concrete is determined by parameters such as aggregate type, strength grade, and age, and can be obtained according to the "Technical Specification for Testing Building Foundation Piles". When the pile has defects such as diameter reduction, diameter expansion, concrete segregation, mud inclusion, or broken piles, the wave velocity value of the defect area will deviate significantly from the reference value of intact concrete. By using a preset deviation threshold as a quantifiable discrimination criterion, qualitative defect identification can be transformed into quantitative spatial area identification. Compared to purely qualitative identification methods that do not introduce reference values and deviation thresholds, the quantitative defect identification of this invention provides a quantifiable criterion for defect judgment. The defect identification results are repeatable and comparable, and the severity of defects can be graded according to different deviation threshold values. This provides a quantitative judgment basis that meets the requirements of current standards for the issuance of quality and safety inspection and testing reports for bridge infrastructure.
[0028] Furthermore, this invention establishes two core subsystems: a probe unit and a surface data processing workstation, and specifies the following three key technical configurations: First, a coupling airbag is installed on the circumferential direction of the probe unit's outer wall. The probe integrates three sensors: a three-component accelerometer, an inclinometer, and a locator. These three sensors share the same time reference and are configured to synchronously output elastic wave signals, three-dimensional spatial attitude parameters, and spatial position information. Second, after the probe is lowered to each observation node, the processor controls the coupling airbag to inflate until it adheres to the borehole wall, then initiates the iteration step corresponding to that observation node. Within the iteration step, the processor reads the attitude parameters output by the inclinometer, calculates the attitude confidence factor, and sends it to the downstream calculation stage. Third, the processor does not accept data from the next observation node until it has completed the calculation of the model update for this iteration. The coupling airbag uses a high-polymer flexible elastic damping material to replace the traditional mechanical pushing device. It can adapt to drilling with different hole diameters under continuously adjustable inflation pressure from 0.1 to 0.6 MPa, achieving uniform fit to the undulations of the hole wall. The three sensors share the same time reference (error less than one microsecond) to ensure a strict temporal correspondence between the elastic wave signal and the probe attitude and spatial position at the same moment, avoiding time misalignment between attitude correction and spatial positioning. The strict control timing of inflation before iteration and no new data being accepted before the iteration is completed ensures that the inversion iteration process will not degrade due to data aliasing. Compared to existing mechanical pushing devices, this coupling airbag, as a structural innovation of the borehole-deployed elastic wave detection probe, can improve the signal-to-noise ratio of signal acquisition by about 6 to 10 dB and eliminate the mode shape differences caused by uneven adhesion. Compared to sensor combinations without a unified time reference, the time sequence correspondence of the synchronous output of the three sensors can ensure the effectiveness of attitude correction and avoid false attitude errors introduced by time misalignment. Compared to existing devices without strict time constraints, this strict time constraint ensures the data consistency of the inversion process of multiple observation nodes, significantly improving the stability of the inversion results and providing hardware guarantee with high signal-to-noise ratio and high time accuracy for regular quality and safety inspection and testing services of bridge infrastructure.
[0029] Furthermore, the real-time data link between the probe unit and the data processing workstation of the present invention is implemented using optical fiber communication. Simultaneously, a storage module is further provided within the probe unit. This storage module is configured to temporarily store the elastic wave signal, the three-dimensional spatial attitude parameters, and the spatial position information when the real-time data link transmission is obstructed. After the real-time data link is restored, the temporarily stored data is transmitted to the data processing workstation via optical fiber communication. Optical fiber communication has advantages over traditional wired cable communication in three aspects: large transmission bandwidth, strong anti-electromagnetic interference capability, and long transmission distance. It is particularly suitable for complex electromagnetic environments such as long signal transmission distances during borehole lowering operations, potential conductivity of fluids within the borehole, and potential electromagnetic interference sources near the surface. The storage module allows the probe to temporarily store the collected data locally when the optical fiber link is obstructed due to mechanical vibration, loose connectors, temporary faults, etc., and then transmit it back after the link is restored, ensuring that the data acquisition process is not interrupted or lost due to communication failures. Compared to existing devices that only use wired cable communication, fiber optic communication significantly improves the signal-to-noise ratio and reliability of data transmission. The distortion of elastic wave signals during transmission is reduced from about 1.5% in wired cable methods to less than 0.1% in fiber optic methods. Compared to existing devices without a temporary storage mechanism, the storage module's temporary storage and delayed transmission capabilities provide link fault tolerance for on-site data acquisition, avoiding the failure of the entire test due to temporary communication failures. This increases the data integrity of the entire test from about 92% to nearly 100%, providing highly reliable data transmission guarantees for quality, safety, and non-destructive testing services for bridge infrastructure in complex on-site environments.
[0030] Furthermore, the electronic device form of this invention expands the industrial application scenarios of this method, including cloud inversion service (centralized inversion of multi-node pile foundation test data), edge computing node inversion (real-time on-site inversion), local workstation inversion (offline in-depth analysis by professional testing personnel), and other deployment methods, providing a flexible technical carrier for the industrial promotion of non-destructive testing services for bridge infrastructure quality and safety.
[0031] This invention has achieved the following systematic beneficial effects at the engineering application level: In terms of applicability, this invention lowers the probe to the pile side borehole for lateral excitation and three-component signal acquisition, fundamentally avoiding the dependence of existing acoustic transmission methods on pre-embedded acoustic tubes, the requirements of low strain methods on the free end of the pile top, and the cost requirements of cross-hole elastic wave CT methods for drilling on both sides. This makes the method applicable to pile foundation testing under various complex conditions such as existing bridges, completed bridges, and pile tops permanently covered by pile caps and piers. It is also applicable to various professional testing service scenarios such as routine quality and safety inspection and testing services for bridge infrastructure, periodic inspection and testing of bridge engineering, and special inspection and testing of bridge pile foundations under complex service conditions.
[0032] In terms of accuracy, this invention achieves a synergistic improvement in accuracy compared to conventional full-waveform inversion without attitude modulation, by reducing the number of inversion iterations by about 60%, reducing the defect location error by about 65%, and increasing the inversion success rate by about 17 percentage points under exemplary working conditions through adaptive inversion driven by attitude confidence factors and dual wave velocity identification based on transverse and longitudinal wave decoupling. This provides a quantitative testing basis with engineering reliability for the quality level assessment, safety status evaluation, and service performance testing of bridge infrastructure.
[0033] In terms of robustness, this invention employs multiple safeguard mechanisms, including an additional attenuation term for the attitude change rate, an attitude confidence factor throughout the inversion iteration, and strict timing constraints on the processor. These mechanisms ensure stable convergence and reliable defect identification even under complex conditions such as a 10 dB signal-to-noise ratio, sudden attitude changes, and temporary communication failures. This meets the accuracy requirements for quality and safety inspection and testing of bridge infrastructure under complex on-site environments and various external interference conditions.
[0034] In terms of engineering convenience, the coupling airbag of this invention allows the device to adapt to drilling with different hole diameters, eliminating the need for customized mechanical pushing devices for different hole diameters. The testing time for a single pile foundation is reduced from several hours in the traditional cross-hole method to less than 30 minutes in this method. The fiber optic communication method and the temporary storage mechanism of the storage module provide link fault tolerance for on-site data acquisition, significantly reducing the on-site operation cost and time consumption for regular quality and safety inspection and testing services for bridge infrastructure, and improving the service efficiency and service capabilities of relevant inspection and testing institutions.
[0035] In summary, the overall technical solution defined by this invention has unexpected technical effects that are significantly superior to existing technologies in engineering applications. It can effectively solve the technical problems of non-destructive testing of existing bridge foundation piles, and provides a high-precision, high-reliability, and engineering-feasible core technical solution for the construction and upgrading of the non-destructive testing and inspection service system for bridge infrastructure quality and safety. It has important engineering value for ensuring the safe operation of bridge infrastructure and maintaining the sustainable service of transportation infrastructure. Attached Figure Description
[0036] Figure 1 A schematic diagram of the overall structure of the borehole elastic wave pile foundation testing device; Among them, 1-drill hole, 2-pillar, 3-ground surface, 4-coupling airbag, 5-pile foundation, 6-steel wire rope.
[0037] Figure 2 This is a schematic cross-sectional view of the internal structure of the probe unit; Figure 3 This is a schematic diagram of hardware data flow and processor control timing. Figure 4The main flowchart for tilt angle adaptive full waveform inversion iteration is shown below. Figure 5 The attitude confidence factor Gaussian kernel mapping characteristic curve; Figure 6 A schematic diagram comparing the convergence speed of the full waveform inversion iteration under conditions of attitude confidence factor modulation and non-attitude confidence factor modulation; Figure 7 A schematic diagram of a network architecture for learning pose confidence factor mapping in a neural network; Figure 8 This is a schematic diagram of a self-supervised attitude-waveform joint correction network architecture; Figure 9 This is a schematic diagram showing the inversion results of transverse and longitudinal wave velocities and the identification of pile defects. Detailed Implementation
[0038] The present invention will be further described in detail below with reference to the accompanying drawings and specific embodiments to facilitate a clear understanding of the present invention, but these descriptions do not constitute a limitation on the present invention.
[0039] Example 1 This embodiment of the device illustrates a typical engineering scenario, used to inspect the reinforced concrete bored piles of an existing cross-river highway bridge. The pile diameter is 1.2 meters, the designed pile length is 35 meters, the pile top is permanently covered by the abutment and piers, and the bridge is in normal traffic operation. The device consists of two main parts: a surface system and a lowering system. (See attached diagram.) Figure 1 The surface system includes a winch, a data processing workstation, and a lateral vibration source hammer; the lowering system is a probe unit lowered by the winch into the lateral borehole of the pile foundation, which is about 1.5 meters from the pile side surface, with a diameter of 110 mm and a depth of 38 meters.
[0040] The winch is an electrically controlled wire rope winch. The wire rope is connected to the probe unit via a load-bearing retainer. Two parallel channels, an optical fiber communication cable and a gas delivery pipeline, are integrated axially within the wire rope. The optical fiber communication cable uses single-mode 9 / 125-micron quartz optical fiber with a polyurethane pressure-resistant sheath. It has a maximum lowering depth of 100 meters and a unidirectional transmission bandwidth of no less than 1Gbps, capable of simultaneously carrying multiple data streams, including elastic wave signals, attitude parameters, spatial position information, and processor downlink control commands. The gas delivery pipeline is a 6mm inner diameter pressure-resistant nylon tube with a pressure rating of 2 MPa, connecting a surface compressed air source to the internal micro-pump gas source interface of the probe unit. The micro-pump acts as a precise pressure regulating unit, adjusting the low-pressure gas source (approximately 0.2 MPa) delivered from the surface compressed air source via the gas delivery pipeline to the required working pressure (0.2 to 0.6 MPa continuously adjustable) for the coupling airbag, and achieving closed-loop pressure regulation through PWM duty cycle control.
[0041] The probe unit adopts a sealed pressure chamber structure, see [link / reference] Figure 2 Its outer shell is made of 304 stainless steel, with an outer diameter of 95 mm and a length of 1.2 meters, and a pressure resistance rating of no less than 1.5 MPa. The interior of the cabin integrates, along the axial direction, a main control module, a storage module, a preprocessing module, a battery, a three-component accelerometer, a high-precision inclinometer, a positioning device, and a miniature air pump. These modules are interconnected via an internal CAN bus and a time synchronization signal line, ensuring that all sensors share the same time reference, and that the timestamp consistency error of the three types of data does not exceed one microsecond.
[0042] The three-component accelerometer is a capacitive triaxial accelerometer fabricated using MEMS technology. It has a measurement range of ±10g, a frequency band of 1 Hz to 5 kHz, and a noise density better than 25 microg / Hz. Its three sensing axes are orthogonally mounted in the X, Y, and Z axes of the probe's coordinate system, used for real-time synchronous acquisition of elastic wave signals during pile-side seismic source excitation. The elastic wave signal (i.e., acoustic signal) acquired by the three-component accelerometer undergoes power frequency interference removal, anti-aliasing filtering, and 24-bit analog-to-digital conversion by the preprocessing module. It is then transmitted to the main control module via the internal CAN bus. The main control module encapsulates the data and uploads it in real-time to the data processing workstation via the fiber optic communication cable.
[0043] The high-precision inclinometer is also a MEMS integrated device, incorporating a dual-mode fusion architecture of a three-axis gyroscope and a three-axis accelerometer. Its inclinometer measurement accuracy is better than 0.05 degrees, azimuth measurement accuracy is better than 0.2 degrees, and sampling rate is 1 kHz. It outputs three-dimensional Euler angle parameters, namely inclinometer component, azimuth component, and rotation component. The positioning device calculates the current lowering depth of the probe based on a fiber optic pressure sensor and a cable tension sensor, with a depth measurement accuracy better than 0.05 meters. It outputs the spatial position x of the probe. The output data from both the inclinometer and the positioning device are collected by an internal CAN bus to the main control module, sharing a unified time reference with the elastic wave signal output from the preprocessing module, ensuring a strict time-series correspondence between the elastic wave signal and the probe's attitude parameters and spatial position information at the same moment.
[0044] The main control module employs a 1.4 GHz quad-core processor based on the ARM Cortex-A53 architecture, running a real-time operating system and handling tasks such as data acquisition scheduling, attitude parameter calculation, air pump control signal generation, and bidirectional communication with the surface workstation. A bidirectional real-time control link is maintained between the main control module and the data processing workstation processor via a fiber optic communication cable. The main control module receives high-level control commands from the data processing workstation processor via the fiber optic communication cable, including inflation commands, acquisition trigger commands, and release permission commands, and converts them into low-level control signals acting on the three-component accelerometer, the inclinometer, the positioning device, and the miniature air pump. Simultaneously, the main control module reports probe status information to the data processing workstation in real time, including airbag pressure confirmation, acquisition completion flags, and data link status. The storage module uses 256 gigabytes of industrial-grade solid-state storage for local temporary storage of all raw data from this construction site. When the transmission of the optical fiber communication cable is obstructed, the main control module temporarily stores the elastic wave signal, the three-dimensional spatial attitude parameters, and the spatial position information in the storage module. After the transmission of the optical fiber communication cable is restored, the main control module reads the temporarily stored data from the storage module and transmits it back to the data processing workstation in sequence via the optical fiber communication cable. The preprocessing module is an FPGA-based hardware signal conditioning unit, responsible for power frequency interference removal, anti-aliasing filtering, and 24-bit analog-to-digital conversion of the three-component signal.
[0045] The battery powers all modules and uses a high-capacity lithium-ion battery pack with a nominal capacity of 50 Ah, a nominal voltage of 14.8 V, and a continuous working time of no less than 24 hours. The coupling airbag is located circumferentially on the outer wall of the probe unit. The coupling airbag is seamlessly stitched with a flexible, elastic damping polymer material (such as a composite braided layer of neoprene rubber and aramid fiber). In its deflated state, its outer diameter is 95 mm, flush with the outer diameter of the probe body. When inflated, it can expand radially to a continuously adjustable range of 130 mm to 200 mm to adapt to different apertures and aperture wall undulations, with a maximum working pressure of 0.6 MPa.
[0046] The micro air pump is a DC brushless piston micro air pump with a maximum output pressure of 0.8 MPa and a maximum flow rate of 5 liters per minute. As a precise pressure regulating unit, it pressurizes the compressed air source from the ground through the gas delivery pipeline to the working pressure required by the coupling airbag, and achieves closed-loop pressure regulation through PWM duty cycle control. An inflation passage is provided between the micro air pump and the coupling airbag. The micro air pump is configured to receive the underlying inflation control signal from the main control module and inflate the coupling airbag until the coupling airbag is in contact with the borehole wall. The underlying inflation control signal of the main control module is triggered and generated by the inflation command sent by the processor of the data processing workstation through the fiber optic communication cable.
[0047] The data processing workstation is a high-performance computing server equipped with a multi-core CPU and a graphics processing unit. The processor is an Intel Xeon Gold 6248R (24 cores, 48 threads, 3.0 GHz), with 256 gigabytes of DDR4 ECC memory and an NVIDIA RTX A6000 graphics processing unit (48 gigabytes of video memory), used to handle the parallel computing tasks of the full waveform inversion iterative process. A real-time data link is established between the workstation and the probe unit via a fiber optic transceiver located on the ground and another fiber optic transceiver located inside the probe unit, connected by the aforementioned fiber optic communication cable. The fiber optic transceivers at both ends perform bidirectional conversion between electrical and optical signals.
[0048] The processor of the workstation strictly follows the hardware control timing sequence when executing the detection method, see [link / reference]. Figure 3 Specifically, after the winch lowers the probe unit to the current observation node (typically, the observation node depth interval is 0.5 meters; for the 38-meter borehole depth in this embodiment, a total of 76 observation nodes are set), the processor first issues an inflation command to the micro air pump, which inflates the coupling airbag until it is tightly attached to the borehole wall (typical inflation time is 3 to 5 seconds). After the airbag attachment status is confirmed, the surface excitation hammer excites an elastic wave source at the top position of the pile side (typical excitation energy is about 200 joules, frequency band covers 100 Hz to 5 kHz). The three-component accelerometer continuously records the elastic wave signal for 200 milliseconds at a sampling rate of 10 kHz, the inclinometer synchronously outputs the three-dimensional spatial attitude parameters of the current observation node, and the positioning instrument outputs the current depth. These three types of data are encapsulated by the main control module and uploaded to the data processing workstation in real time through the fiber optic communication cable. After receiving complete data from the observation node, the workstation's processor initiates the corresponding full-waveform inversion iteration step. Only after all calculations in this iteration step are completed does it issue a "allow lowering to the next observation node" command to the surface winch controller. The workstation controls the inflation and deflation of the coupling airbag by controlling a miniature air pump and controls the position of the probe unit by controlling the length of the wire rope on the winch. This strict timing ensures that the inversion iteration process will not introduce data aliasing errors due to the previous observation node not being fully processed.
[0049] Example 2 This method embodiment inherits the hardware configuration of the aforementioned device embodiment and conducts a complete inspection process for the foundation piles of the cross-river highway bridge. The entire inspection process is divided into four main stages: on-site data acquisition, data preprocessing and coordinate correction, full waveform inversion iteration, defect identification, and report output. See [link to documentation]. Figure 4During the inversion iteration stage, preferred implementation methods based on neural networks and reinforcement learning are interspersed after the attitude confidence factor calculation, adaptive step size factor calculation, and gradient direction calculation stages, respectively, to fully demonstrate the feasibility and technical effectiveness of the present invention.
[0050] Step S1: Acquisition of elastic wave signals from multiple observation nodes.
[0051] Following the hardware deployment and airbag coupling control process described in the aforementioned device embodiment, the on-site data acquisition team set up 76 observation nodes at 0.5-meter intervals along the depth of the pile foundation in the pile-side borehole, collecting data section by section from the bottom of the hole (38 meters deep). At each observation node, five actions were sequentially performed: probe airbag inflation, ground surface vibration, three-component signal recording, synchronous acquisition of attitude parameters, and spatial position information recording. The data acquisition for all 76 observation nodes took approximately 2.5 hours.
[0052] At each observation node, the three-component accelerometer acquires elastic wave signals in real time. These signals are then uploaded to the data processing workstation via the aforementioned fiber optic communication link, serving as the observation data for subsequent full waveform inversion. In this system, the subscript 's' represents the source number, the subscript 'r' represents the detector number, the superscript 'obs' represents the observation data, and 't' represents time. Simultaneously, the inclinometer outputs the three-dimensional spatial attitude parameters of the probe at the current observation node, namely the tilt angle θ(t), azimuth angle φ(t), and rotation angle ψ(t), while the locator outputs the spatial position x of the probe at the current observation node. All acquired elastic wave signals, attitude parameters, and spatial position information carry a timestamp based on the same time base, facilitating strict time alignment in subsequent algorithmic steps.
[0053] Step S2, geometric preprocessing correction based on attitude parameters.
[0054] Since the actual attitude of the probe at different observation nodes may be slightly deflected relative to the ideal vertical attitude due to factors such as irregular borehole walls and asymmetrical airbag expansion (typical deflection angle is between ±3 degrees and ±8 degrees), the signals collected by the three-component accelerometer are actually signals in the probe body coordinate system. They need to be transformed to the absolute geographic reference coordinate system in order to correctly participate in the subsequent full waveform inversion based on absolute geometry.
[0055] Specifically, the processor of the workstation constructs a three-dimensional spatial rotation matrix based on the attitude parameters acquired at each observation node. This rotation matrix is obtained by multiplying the rotation matrices around the Z-axis, the Y-axis, and the X-axis in the order of Euler angles ZYX. Then, the processor performs an inverse transformation of the spatial coordinate system on the original three-component signal of the observation node to obtain a corrected three-component signal, which is then used as the observed waveform in the subsequent objective function construction.
[0056] Preferably, based on the geometric preprocessing correction, the processor of the workstation can also execute an attitude-frequency coupled domain S-transform correction method as a preferred embodiment to address the complex working condition where attitude changes during probe descent exhibit differentiated distributions across different frequency bands. In the low-frequency band (typically less than 30 Hz), attitude changes are dominated by slow drift, corresponding to changes in winch descent speed and accumulated contact friction with the borehole wall; in the mid-to-high frequency band (typically 30 Hz to 200 Hz), attitude changes exhibit mixed characteristics, including airbag inflation-wall contact micro-vibration and cable coupling vibration; in the high-frequency band (typically greater than 200 Hz), attitude changes are dominated by high-frequency vibration, corresponding to probe housing resonance and sensor's own mechanical noise. The single rotation matrix used in the geometric preprocessing correction cannot handle the attitude deviations in the above three frequency bands separately. This embodiment decomposes the signal into the time-frequency domain using S-transform, constructs a frequency-dependent rotation matrix for each time-frequency atom, and performs rotation correction, thereby achieving fine correction in the attitude-frequency coupled domain.
[0057] The processor performs an S-transform on the three-component original signal x(t) at each observation node to obtain the time-frequency atom matrix S(t,f). The specific expression of the S-transform is as follows: Where: S(t,f) represents the time-frequency atom value of signal x(t) at time t and frequency f, which is a complex number; x(t) represents the original signal (referring to a component of the three-component acceleration signal, obtained by performing three S-transforms on the x, y, and z components respectively to obtain three time-frequency atom matrices); τ represents the integration variable, with the same dimension as time t; |f| / √(2π)·exp(-(τ-t) 2 f 2 / 2) represents the frequency-dependent Gaussian window function, whose window width adapts to the frequency. At low frequencies, the window width is large and the frequency resolution is high, while at high frequencies, the window width is small and the time resolution is high. exp(-j2πfτ) represents the complex exponential kernel function of the Fourier transform, which is used to extract the phase information at frequency f. j represents the imaginary unit.
[0058] The processor performs frequency band division on the obtained time-frequency atom matrix, dividing it into three frequency bands: the first band is 0.5 Hz to 30 Hz (low-frequency attitude drift dominant band), the second band is 30 Hz to 200 Hz (mixed band), and the third band is 200 Hz to 1000 Hz (high-frequency vibration dominant band). A cosine square smoothing transition function is used at the frequency band boundaries (30 Hz and 200 Hz) to avoid spectral breakage. Specifically, the smoothing transition function takes the following form: at the boundary frequency f... b Within the ±5 Hz transition range, the weight of the low-frequency side is taken as cos 2 (π(ff b +5) / 20), the weights on the high-frequency side are taken as sin 2 (π(ff b +5) / 20), the sum of the two is always equal to 1.
[0059] The processor constructs a frequency-dependent rotation matrix R(α,t,f) at each time-frequency atom (t,f). The attitude parameter α(t,f) of the rotation matrix is obtained by fusing two branches: the first branch faces the low-frequency band (less than 30 Hz), and the attitude parameter is taken from the attitude low-pass component obtained by low-pass filtering the output of the inclinometer with a cutoff frequency of 10 Hz. This component mainly reflects the attitude drift during the winch descent process. The second branch faces the mid-to-high frequency band (greater than 30 Hz), and the attitude parameter is taken from the instantaneous polarization principal direction obtained by three-component polarization principal component analysis (i.e., performing eigenvalue decomposition on the three-component complex covariance matrix at this time-frequency atom). This direction reflects the instantaneous attitude under high-frequency vibration conditions. In the mixed band from 30 Hz to 200 Hz, the attitude parameters of the two branches are weighted and fused according to the aforementioned cosine square smoothing transition function.
[0060] The processor performs a frequency-dependent rotation correction operation on each time-frequency atom matrix to obtain the corrected time-frequency atom matrix: Where X(t,f) represents the complex value vector of the original three-component signal at the time-frequency atom (t,f). R(α,t,f) represents the column vector composed of the complex values of the three component S-transform results at the time-frequency atom; R(α,t,f) represents the three-dimensional spatial rotation matrix constructed at the time-frequency atom (t,f) based on the fused attitude parameters of the two branches; X'(t,f) represents the time-frequency atom vector after rotation correction; the superscript T indicates matrix transpose. The processor performs the above rotation correction operation on each time-frequency atom.
[0061] The processor performs an inverse S-transform on all corrected time-frequency atom matrices to obtain a frequency-dependent decoupled three-component signal: Where: u'(t) represents the frequency-dependent decoupled three-component signal, which is used as the observed waveform in the subsequent objective function construction; S -1 {·} represents the inverse S-transform operator, which specifically involves first integrating X'(t,f) along the frequency dimension to obtain the time-domain signal, and then obtaining an approximation of the original time-domain signal through phase reconstruction and amplitude correction. In engineering practice, the inverse S-transform is numerically implemented using a discrete S-transform based on the fast Fourier transform. A single time-frequency analysis-correction-inverse transform cycle takes approximately 30 milliseconds, which is negligible compared to the overhead of a single step in the inversion iteration (approximately 8 seconds).
[0062] This embodiment, as a frequency-dependent enhancement scheme for geometric preprocessing correction, can form a collaborative working mode with geometric preprocessing: when the processor detects that the attitude distribution of the current observation node satisfies the Gaussian assumption (specifically, the attitude parameter sequence of the node has a p-value greater than 0.05 after the Shapiro-Wilk normality test) and the main energy of the attitude change spectrum is concentrated in a single frequency band (specifically, the energy concentration index of the attitude change rate spectrum exceeds 0.85), geometric preprocessing correction can meet the accuracy requirements; when the processor detects that the attitude distribution deviates significantly from the Gaussian assumption or the attitude change spectrum spans multiple frequency bands, it automatically switches to the attitude-frequency coupled domain S-transform correction described in this extended embodiment. In engineering implementation, the processor first performs the above automatic determination at each observation node, and then selects the corresponding correction method to achieve intelligent switching between the two methods.
[0063] The control experiment data for this embodiment are as follows: In complex conditions where low-frequency attitude drift (typical frequency less than 5 Hz, drift amplitude 0.3 degrees) and high-frequency vibration (typical frequency greater than 200 Hz, vibration amplitude 0.05 degrees) coexist, the single rotation matrix used in the geometric preprocessing correction described above cannot handle the attitude deviations of the two frequency bands separately. This results in the corrected three-component polarization principal axis deviating from the source-detection direction by a maximum of approximately 1.2 degrees, and the velocity inversion deviation being approximately 6.0%. This embodiment performs frequency-dependent rotation correction on each time-frequency atom individually in the time-frequency domain using S-transform, compressing the deviation of the three-component polarization principal axis to within 0.3 degrees and reducing the velocity inversion deviation to 1.0%. In terms of inversion success rate, this embodiment improves the inversion success rate to over 96%. The above comparative experiments fully demonstrate that this embodiment, through fine correction in the attitude-frequency coupling domain, achieves a significant improvement in the ability to handle frequency-dependent attitude deviations, providing higher-quality observation data for subsequent full-waveform inversion.
[0064] Step S3: Forward modeling and objective function construction based on the isotropic elastic wave equation.
[0065] The workstation establishes an initial three-dimensional underground model of the pile foundation and the surrounding medium. This underground model is characterized by model parameters m. In this embodiment, the model parameters are parameterized as three spatial distribution fields: the first Lamé coefficient λ, the second Lamé coefficient μ, and the density ρ of the isotropic elastic medium. The initial model is constructed based on the pile diameter, pile length, concrete strength grade (C35) provided in the pile foundation design documents, and the layered information of the surrounding soil layers provided in the geological survey report. The wave velocity of the pile body in the initial model is taken from the design value, with a longitudinal wave velocity of approximately 4000 m / s and a transverse wave velocity of approximately 2300 m / s. The wave velocity of the surrounding soil layers is taken according to the layered values in the geological survey report.
[0066] Forward modeling is based on the isotropic elastic wave equation, the specific expression of which is: Where: ρ(x) represents the medium density at spatial coordinate x, u(x,t) represents the displacement vector, σ(x,t) represents the stress tensor at spatial coordinate x at time t, f(x,t) represents the volumetric force vector per unit volume at spatial coordinate x at time t, including the excitation force from the external source, and ∇· represents the divergence operator; x represents the spatial coordinate vector, and t represents time. The equations are solved discretized in the time-space domain using the finite difference method, with a time step of 1 microsecond and a spatial grid step of 5 centimeters. The forward modeling calculation takes approximately 8 seconds per source excitation. The wavefield obtained from the forward modeling simulation is denoted as the forward modeling waveform, which serves as the prediction data in the subsequent construction of the objective function.
[0067] The processor of the workstation uses the waveform difference between the observed data and the forward simulation data as the objective function, and the specific expression of the objective function J(m) is: Where: J(m) represents the objective function; m represents the model parameters; N s N represents the total number of earthquake focal points; r s represents the total number of detectors; r represents the number of seismic sources; s represents the number of detectors. Let represent the observed waveform of the elastic wave signal excited by the s-th source, received by the r-th detector at time t. This represents the forward-modeled waveform of the elastic wave signal generated by the s-th source at time t, obtained by forward modeling the isotropic elastic wave equation under model parameters m; T represents the total duration of waveform recording after a single source excitation; t represents time. 2 L represents the vector 2 The square of the norm is the sum of the squares of each component of the three-component signal.
[0068] Step S4: Tilt angle adaptive full waveform inversion iteration.
[0069] This stage is the core algorithm of the invention. The processor of the workstation uses the sensor attitude information measured by the tiltmeter to construct an adaptive update operator, dynamically adjusting the model update direction and step size during the inversion process. Specifically, each inversion iteration sequentially executes five sub-steps: attitude confidence factor calculation, gradient calculation, adaptive step size factor calculation, model update direction construction, model parameter update, and convergence criterion verification. The preferred implementation methods for each sub-step and the preferred implementation method integrating neural networks and reinforcement learning are described below.
[0070] Step 4.1. Calculate the attitude confidence factor.
[0071] The processor extracts tilt, azimuth, and rotation components from the attitude parameters of the observation node corresponding to the current iteration step, and constructs a tilt feature vector: Where: α(t) represents the tilt eigenvector; θ(t) represents the tilt component; φ(t) represents the azimuth component; ψ(t) represents the rotation component; and the superscript T indicates matrix transpose. Then, the processor, centered on the pre-calibrated ideal tilt-free probe attitude α0 (in this embodiment, α0 corresponds to the ideal attitude when the tilt angle is 0, i.e., the tilt component, azimuth component, and rotation component are all 0), uses a Gaussian kernel function to nonlinearly map the difference between the tilt eigenvector and the ideal calibration attitude into an attitude confidence factor. The specific expression for the attitude confidence factor η(t) is: Where: η(t) represents the attitude confidence factor; α(t) represents the aforementioned tilt angle feature vector; α0 represents the calibration attitude, corresponding to the ideal attitude when the tilt angle is 0; This represents the tilt confidence bandwidth parameter, used to control the sensitivity of attitude deviation to confidence attenuation. In this embodiment, it is pre-determined based on approximately three times the laboratory calibration noise standard deviation of this model of tiltmeter, which is approximately 0.05 radians. It equals 0.15 radians. The characteristic curve of the mapping function is as follows: Figure 5 As shown, the attitude confidence factor η(t) takes the maximum value of 1 when the attitude deviation is 0 (corresponding to the ideal vertical attitude and the observation data is completely reliable), decreases monotonically with the attitude deviation, and approaches 0 when the deviation exceeds 3 times the bandwidth parameter (corresponding to a serious deviation from the ideal attitude and the observation data is almost unreliable).
[0072] In a preferred embodiment, the Gaussian kernel function used in the aforementioned attitude confidence factor calculation step can be replaced by a nonlinear mapping of a pre-trained neural network, i.e., the mapping is performed by the pre-trained neural network f. θ Output attitude confidence factor η(t), see Figure 7This preferred embodiment only replaces the specific calculation method of the attitude confidence factor described above; all other steps (gradient calculation, step size factor calculation, model update direction construction, model parameter update, convergence criterion verification) are completely consistent with the previous embodiment.
[0073] In terms of network architecture, f θ A lightweight architecture based on a Transformer encoder is adopted. The input layer contains four sets of features: the tilt feature vector α(t) of the current iteration (3D), the calibrated pose α0 (3D), the normalized iteration progress (1D), and the pose sequence of the five neighboring observation nodes before the current observation node (forming a 15-dimensional sequence feature to capture the temporal correlation of pose changes). These features are first position-encoded and then fed into the Transformer encoder backbone. The backbone consists of two stacked standard Transformer encoder layers, each containing a 4-head self-attention module (embedding dimension equal to 64, key dimension equal to 16) and a feedforward network module (hidden dimension 128, activation function GELU). Residual connections and layer normalization are set between layers to stabilize training. The Transformer encoder output is fed into the output fully connected layer (64 to 32 to 1, the last layer uses the Sigmoid activation function to limit the output to the range of 0 to 1) after global average pooling, and finally outputs a scalar pose confidence factor η(t), whose value range is consistent with η(t) in the previous implementation.
[0074] The training data construction employs a two-stage strategy combining synthetic pile foundation inversion datasets with fine-tuning using engineering measured data. The first stage uses synthetic data for pre-training, constructing a training sample containing approximately 5000 different pile foundation models (covering typical pile foundation conditions such as intact piles, reduced-diameter piles, expanded-diameter piles, broken piles, concrete segregation piles, and mud-filled piles; pile lengths are randomly distributed between 15 and 50 meters, pile diameters are randomly distributed between 0.6 and 2.0 meters, and the surrounding soil structure is randomly generated according to the actual geological report of the project). Elastic wave signals are simulated at 76 observation nodes on each pile foundation model, and random perturbations are injected into the attitude parameters of each observation node to simulate attitude deviations in actual engineering (the perturbation amplitude follows a normal distribution with a mean of 0 and a standard deviation of 0.1 radians). The second stage uses engineering measured data for fine-tuning. Historical real pile foundation inversion data whose pile integrity status has been confirmed through destructive testing are used as fine-tuning samples to further fine-tune the network parameters with a small learning rate to adapt to specific engineering geological conditions.
[0075] The loss function design adopts a meta-learning paradigm: the squared L2 norm of the difference between the converged model parameters and the true model parameters, plus the weighting coefficients multiplied by the accumulated residuals of the objective function J(m) across all iterations during the inversion iteration, is used as the outer meta-objective. This outer meta-objective calculates the gradients of the network parameters and updates them through backpropagation based on implicit differentiation. Since the backpropagation computation for a complete inversion iteration is computationally intensive, a truncated BPTT strategy is adopted in practice, retaining only the gradient history of the most recent 20 steps to balance training efficiency and accuracy. The training configuration uses the Adam optimizer with an initial learning rate of 1×10⁻⁶. -4 With a batch size of 32 and 200 training rounds, a single complete training session takes approximately 36 hours on the aforementioned hardware configuration.
[0076] Regarding the inference process, the network parameters are frozen after training. During inference, each inversion iteration requires only one forward computation (taking approximately 2 milliseconds), which is negligible compared to the overhead of a single adjoint state gradient calculation (approximately 8 seconds). The processor directly calls the pre-trained network f during the attitude confidence factor calculation stage of each inversion iteration. θ Output the attitude confidence factor η(t), and the rest of the steps are consistent with the previous implementation method.
[0077] The comparative experimental data for this preferred embodiment are as follows: Among the 76 observation nodes of the aforementioned cross-river highway bridge pile foundation, the attitude deviation at approximately 18 observation nodes exhibited a non-Gaussian distribution characteristic (in actual engineering, irregular borehole walls caused slight jamming-slippage type attitude abrupt changes in the probe, resulting in a bimodal attitude distribution that deviated significantly from the Gaussian kernel assumption). Under this condition, the neural network mapping of this preferred embodiment, by learning the actual attitude distribution characteristics in historical data, can adaptively output a more reasonable attitude confidence factor for this type of bimodal attitude distribution, further reducing the number of inversion iterations from 80 in the previous embodiment to 55, the defect location error from 0.12 meters to 0.07 meters, and expanding the noise robustness boundary from a signal-to-noise ratio of 10 dB to 6 dB.
[0078] Step 4.2, calculate the gradient.
[0079] The processor calculates the objective function J(m) with respect to model parameters m under the current model parameters based on the adjoint state method. i The gradient of , wherein the specific expression of the gradient is: in: Let i represent the i-th component of the model parameter m, where i is the index of the model parameter component. Traverse all discrete degrees of freedom of the first Lamé coefficient λ, the second Lamé coefficient μ, and the density ρ on the spatial grid. J(m) represents the aforementioned full waveform inversion objective function, which is consistent with the variable of the same name in step S3; The group represents the forward modeling waveform of the elastic wave signal excited by the s-th source, obtained by forward modeling of the isotropic elastic wave equation under the current model parameters m; This represents the second-order partial derivative of the forward wave field with respect to time t; The accompanying wave field is represented by the difference between the observed waveform and the forward waveform. The equation of the accompanying elastic wave is obtained by injecting it with the source term in reverse time. The subscript i is consistent with the gradient component index i mentioned above, and corresponds to the i-th model parameter component. T represents the total duration of waveform recording after a single source excitation. t represents time.
[0080] It should be noted that the specific integrand form of the gradient expression above is based on the model parameter components. Corresponding to density The situation at that time should be given; For the first Lamé coefficient λ or the second Lamé coefficient μ, the specific integrand form of the adjoint state method gradient involves the product of the spatial derivative (strain tensor) of the forward wavefield and the corresponding component of the adjoint wavefield. Its specific form can be obtained according to the standard derivation of the adjoint state method, and those skilled in the art can implement it as needed. Partial derivatives with respect to the index i of all model parameter components... Together they form a complete gradient vector This gradient vector serves as the input for the gradient direction construction step in step 4.5.
[0081] Step 4.3 Calculate the adaptive step size factor. The processor calculates the tilt angle adaptive step size factor η(k) for the current iteration based on the aforementioned attitude confidence factor η(t), the difference in tilt angle feature vectors between two adjacent iterations, and a preset initial step size. The specific expression for the tilt angle adaptive step size factor is as follows: in: ηk represents the tilt angle adaptive step size factor for the current k-th iteration; η0 represents the initial step size, which is 0.05 in this embodiment; N t Indicates the number of time sampling points; t n Represents the nth time sampling point; η(t) n ) represents the attitude confidence factor calculated at the nth time sampling point according to the aforementioned attitude confidence factor expression; γ represents the adaptive adjustment coefficient, which is taken as 5.0 in this embodiment; α (k) α represents the tilt angle eigenvector of the current k-th iteration;(k-1) This represents the tilt angle feature vector from the (k-1)th iteration. When the attitude confidence factor η(t) is close to 1, the tilt angle adaptive step size factor η... (k) The initial step size η0 is kept close to ensure strong search capability; when the attitude confidence factor η(t) approaches 0, the tilt angle adaptive step size factor η (k) The step size approaches zero, adaptively shielding the influence of the observation data from this observation node on the inversion. Furthermore, the exponential decay term of the difference in the tilt eigenvector between adjacent iterations further ensures that the step size is dynamically compressed when attitude parameters fluctuate drastically, thereby improving the stability of the inversion iteration process.
[0082] In another preferred embodiment, the aforementioned tilt angle adaptive step size factor calculation based on the weighted summation of attitude confidence factors can be performed using a policy network π pre-trained through reinforcement learning. φ The preferred implementation replaces the dynamic decision-making process described above. It only replaces the specific calculation method of the tilt angle adaptive step size factor in the previous implementation, while all other steps (attitude confidence factor calculation, gradient calculation, model update direction construction, model parameter update, and convergence criterion verification) remain consistent with the previous implementation.
[0083] The reinforcement learning framework employs the proximal policy optimization (PPO) algorithm. The state space s k Defined as a 7-dimensional vector: the ratio of the current iteration count to the maximum iteration count (1-dimensional), the norm of the current gradient vector (1-dimensional), the residual sequence of the objective function J(m) from the last 3 iterations (3-dimensional), the attitude confidence factor η(t) of the current observation node (1-dimensional), and the attitude change rate (i.e., the norm of the difference in the tilt eigenvectors between two adjacent iterations) (1-dimensional). The action space is a continuous 2-dimensional vector: tilt adaptive step size factor η. (k) The value range [0.001, 0.5] and the momentum term coefficient β (k) The value range is [0.0, 0.99]. The reward function r k It consists of three items: the first reward is the negative of the decrease in the objective function (to encourage faster convergence), the second reward is the coefficient -0.01 multiplied by the gradient oscillation (to encourage stable descent), and the third reward is the final state reward -10 multiplied by the final model parameter bias (to encourage inversion accuracy).
[0084] Policy Network π φ With Value Network V ψAll models employed a 3-layer multilayer perceptron architecture (hidden layers 64 to 64 to the output layer). Training used a PPO-Clip variant with a truncation threshold of 0.2, a generalized dominance estimation parameter of 0.95, a discount factor of 0.99, a single batch trajectory length of 1024, and 16 parallel environments. Training data consisted of approximately 500 synthetic pile inversion episodes (each episode completed one full waveform inversion iteration until convergence or the maximum number of iterations of 200 was reached). Training on the aforementioned hardware configuration took approximately 72 hours.
[0085] Regarding the inference process, the policy network is frozen after training. During inference, the policy network is called at each inversion iteration to output the current tilt angle adaptive step size factor η. (k) With momentum term coefficient β (k) This replaces the tilt angle adaptive step size factor calculation formula in the above embodiments. The tilt angle adaptive step size factor η (k) With momentum term coefficient β (k) It is then fed into the next node's model update direction construction stage.
[0086] The comparative experimental data for this preferred embodiment are as follows: When the inversion iteration enters the saddle point platform (i.e., the objective function gradient is very small but far from the optimal solution, and conventional full waveform inversion is prone to getting trapped in local minima in this state), the reinforcement learning strategy of this preferred embodiment can actively judge the characteristics of the saddle point platform by having a holistic perception of the state space and temporarily increase the step size to jump out of the platform, so that the inversion success rate (i.e. the proportion of reaching the global optimum within 200 iterations) is increased from about 78% in the previous embodiment to about 95%.
[0087] Step 4.4 Constructing the model update direction.
[0088] The processor constructs the adaptive update direction for the current iteration based on the aforementioned tilt angle adaptive step size factor η (k), preprocessing operator P (k), gradient vector, momentum term coefficient β (k), and the model update amount Δm (k-1) from the previous iteration. The specific expression for the adaptive update direction (i.e., the model update amount) is as follows: Where: Δm (k) Indicates the model update amount (adaptive update direction) in the current k-th iteration; k represents the iteration number; η (k) P represents the aforementioned tilt angle adaptive step size factor; (k) The preprocessing operator is represented (in this embodiment, a diagonal pseudo-Hessian operator is used as an approximation to balance the differences in gradient sensitivity in different spatial regions); ∇ m J (k) β represents the gradient vector of the current k-th iteration;(k) Δm represents the momentum term coefficient, which is taken as 0.9 in this embodiment; (k-1) This represents the model update amount in the (k-1)th iteration.
[0089] In another preferred embodiment, the aforementioned geometric preprocessing correction step based on attitude parameters can be implemented using a self-supervised learning-based attitude correction network g. φ Alternative, see Figure 8 The self-supervised learning attitude correction network denoises and smooths the attitude parameters output by the original inclinometer, thereby eliminating the influence of attitude parameter errors caused by factors such as inclinometer measurement noise and airbag expansion vibration on the subsequent coordinate transformation. This preferred embodiment operates in the geometric preprocessing correction stage, and the corrected attitude parameters output thereafter enter the subsequent forward simulation, attitude confidence factor calculation, gradient calculation, adaptive step size factor calculation, and model update direction construction stages.
[0090] In terms of network architecture, g φ An architecture combining Bi-LSTM and fully connected layers is adopted. The input layer receives the original 9-dimensional data sequence from the inclinometer (including the outputs of the three-axis gyroscope, three-axis accelerometer, and three-axis magnetometer, forming a complete 9-DOF attitude measurement data stream) and the initial attitude estimate (the real-time attitude output by the Kalman filter built into the inclinometer). The backbone is a two-layer bidirectional LSTM network with 128 hidden units in each layer, capturing the bidirectional correlation of the attitude sequence in the time dimension. The Bi-LSTM output is mapped to the corrected three-dimensional Euler angles, θ(t), φ(t), and ψ(t), through fully connected layers (256 to 64 to 3), which are consistent with the attitude components directly output by the inclinometer in the previous implementation.
[0091] The design of the self-supervised loss function is a key innovation of this preferred implementation. Specifically, the total loss L... total It consists of a weighted sum of three terms: the first term is the waveform difference loss L. waveform The corrected waveform is obtained by constructing a rotation matrix from the corrected attitude, performing an inverse coordinate transformation on the original three-component signal, obtaining the forward waveform through forward modeling, and calculating the L2-norm squared error between the two, with weight α. loss Equals 1.0; the second term, pose smoothing loss L smooth The L2 norm of the difference in attitude parameters between adjacent iteration steps (used to constrain the physical continuity of attitude changes and avoid abrupt changes in network output), and the weight β. loss Equal to 0.1; Third term polarization principal axis uniformity loss L phys The distance is equal to the cosine of the angle between the principal polarization direction and the source-detector direction at the arrival time of the P-wave head of the corrected three-component signal (used to force the network output to conform to the physical principle of elastic wave polarization), and the weight γ. lossIt equals 0.05. Among the three loss functions, the backpropagation of the waveform difference loss needs to pass through two nontrivial operators: rotation matrix construction and forward simulation. In engineering, an end-to-end implementation based on automatic differentiation is adopted, and the adjoint state method is used to accelerate gradient calculation.
[0092] Regarding training data and procedures, a two-stage strategy of pre-training with synthetic data followed by fine-tuning with measured data is adopted, similar to the preferred implementation of neural network confidence factor mapping. The pre-training stage uses attitude-signal pairing samples from approximately 8000 synthetic pile foundation observation nodes, and the training configuration employs the Adam optimizer with a learning rate of 5×10⁻⁶. -4 Batch size 16, number of training rounds 150.
[0093] In terms of the inference process, after receiving the original 9-dimensional sequence output by the inclinometer, the processor calls the pre-trained network g with index φ to directly output the corrected attitude parameters. These corrected attitude parameters then replace the original attitude parameters in the previous implementation method, and proceed to the subsequent rotation matrix construction and coordinate system inverse transformation process. A single forward calculation takes approximately 5 milliseconds, which has a negligible impact on the overall inversion process.
[0094] The comparative experimental data for this preferred embodiment are as follows: On engineering measured data with 5% inclinometer measurement noise (standard deviation 0.1 degrees), the geometric preprocessing correction of the aforementioned implementation method exhibits significant residual phase crosstalk (the three-component polarization principal axis deviates from the source-detection direction by up to approximately 4.5 degrees), resulting in a velocity inversion deviation of approximately 6% even after subsequent full waveform inversion convergence. The self-supervised attitude-waveform joint correction of this preferred embodiment effectively suppresses the propagation of inclinometer noise through end-to-end learning, compressing the polarization principal axis deviation to within 1.2 degrees and reducing the velocity inversion deviation to 1.8%. The aforementioned three preferred embodiments are used in combination (i.e., first g... φ Correct the attitude parameters, then use f θ Output attitude confidence factor η(t), then use π φ Output tilt angle adaptive step size factor η (k) The inversion accuracy can be further improved by about 25%.
[0095] The workstation's processor can also execute a physical model-based attitude-waveform joint inversion method as another extended implementation. In the previous implementation, the attitude parameters are frozen after one rotation in the geometric preprocessing correction stage, and are not updated during the inversion iteration process. When the attitude measurement itself has significant noise, or the attitude distribution deviates significantly from the Gaussian confidence factor assumption, the attitude error remaining from the first rotation in preprocessing will be transmitted to the inversion process through the observation data, limiting further improvement in inversion accuracy. This preferred implementation elevates the attitude parameter R(α) to an inversion variable, participating in the inversion together with the model parameter m. During the inversion iteration process, the attitude parameters and model parameters are updated collaboratively, achieving secondary correction of the residual attitude error from preprocessing.
[0096] The processor constructs an extended objective function J(m,R), using the attitude parameter R(α) and the model parameter m as variables to be inverted. The specific expression of the extended objective function is as follows: Where: J(m,R) represents the extended objective function; m represents the model parameters (i.e., the spatial distribution field composed of the first Lamé coefficient λ, the second Lamé coefficient μ, and the density ρ); R(α) represents the three-dimensional spatial rotation matrix constructed from the attitude parameter α, where The attitude parameters to be inverted; The meaning is consistent with the corresponding symbol in the aforementioned objective function J(m); The waveform represents the observation (consistent with the observation waveform in the aforementioned implementation, the original three-component observation signal without geometric preprocessing correction). This represents the forward simulation waveform under model parameter m; This represents the predicted waveform obtained by transforming the forward-modeled waveform to the probe body coordinate system using a rotation matrix; Φ(R) represents the prior regularization term of the attitude parameters, specifically in the form of... ,in To calibrate the attitude; The prior weight is represented by R(α), which is set to 0.01 in this embodiment. In this preferred embodiment, the forward modeling waveform is rotated to the probe body coordinate system via R(α) and then compared with the original observed waveform. Thus, the attitude parameter α can be used together with the model parameter m as inversion variables to participate in gradient calculation and iterative update.
[0097] The processor calculates the extended objective function J(m,R) with respect to the model parameter components. With regard to attitude parameter components The gradient. Regarding the model parameter components. The specific expression for the gradient is: It should be added that the accompanying wave field The difference between the observed waveform and the forward-modeled waveform after rotation matrix transformation is used as the basis for solving the elastic wave equation propagating backward in time from the accompanying source. That is, the specific form of the accompanying source is... The accompanying earthquake source in the previous implementation method was... The direct composition differs.
[0098] Regarding attitude parameter components The specific expression for the gradient of (where j is the index of the attitude parameter component, taking the tilt angle, azimuth angle, and rotation angle in sequence) is: in: Represents the rotation matrix R(α) with respect to the attitude parameter components. The partial derivative matrix can be directly obtained from the analytical expression of R(α). Element-wise differentiation yields the result. Since the rotation matrix in this embodiment is a product of ZYX Euler angles,... The analytical form is the product of the derivative of the rotation matrix containing α_j and the other rotation matrices in the aforementioned chain expression; and The meaning is consistent with the corresponding symbol in the aforementioned extended objective function J(m,R); Represents the attitude deviation vector The j-th component; the last term Denotes the prior regularization term Φ(R) for all pairs of ... The partial derivatives. The gradient calculation can be directly implemented using the backpropagation function of automatic differentiation software libraries (such as PyTorch, JAX). In engineering practice, the time for a single attitude gradient calculation is about 0.5 seconds, which is much less than the time for a single model gradient calculation (about 8 seconds).
[0099] The processor can employ two optimization strategies to jointly optimize the extended objective function J(m,R) during the inversion iteration process. In the first optimization strategy (i.e., the alternating optimization strategy, which is the default strategy in this embodiment), the processor executes the following sub-steps in each large inversion iteration: First, fix the attitude parameter α and perform one round of tilt angle adaptive inversion iteration as described in the previous implementation (including five sub-steps: attitude confidence factor calculation, gradient calculation, adaptive step size factor calculation, model update direction construction, model parameter update, and convergence criterion verification); then, fix the updated model parameter m and perform one Gauss-Newton iteration on the attitude parameter α, i.e., according to... Search direction for constructing attitude parameters Then, the attitude parameter update step size is determined according to the Wolfe line search criterion. ,according to Update the attitude parameters; finally, determine the convergence condition for the large iteration (i.e., the relative decrease of the expanded objective function is less than the preset threshold of 1×10). -4 If convergence is not achieved, the next large iteration begins; if convergence is achieved, the inversion terminates. In the second optimization strategy (i.e., the joint descent strategy), the processor constructs a block preconditioner using [m;α] as the unified variable. ,in Preconditioners for model parameters (i.e., preprocessing operators in the previous implementation) The attitude parameters are preconditioned (using the diagonal approximation of their Hessian matrices), and then joint L-BFGS iterations are performed until convergence. This embodiment employs the first optimization strategy to fully reuse the inversion iterative framework of the previous implementation.
[0100] The evolution of this implementation method from geometric preprocessing correction is as follows: The previous implementation method follows a two-stage paradigm of correction followed by inversion, that is, before the inversion begins, the original observation signal is transformed to the absolute geographic reference coordinate system through a rotation during geometric preprocessing. During the inversion iteration, the attitude parameters are frozen, and only the model parameters participate in the iterative update. This extended implementation method follows a unified paradigm of attitude as the inversion variable, that is, during the inversion iteration, the attitude parameters and model parameters are used together as variables to be inverted and participate in the iterative update. The secondary correction of the residual attitude error in the preprocessing is achieved through joint optimization of the extended objective function J(m,R). The two implementation methods can be flexibly selected in engineering: when the attitude measurement noise is small (e.g., the inclinometer accuracy is better than 0.05 degrees) and the attitude distribution satisfies the Gaussian assumption, the previous implementation method can meet the accuracy requirements; when the attitude measurement noise is large (e.g., the inclinometer accuracy is worse than 0.2 degrees) or the attitude distribution deviates significantly from the Gaussian assumption, this extended implementation method can significantly improve the inversion accuracy.
[0101] The comparative experimental data for this preferred embodiment are as follows: Under conditions where the attitude measurement noise is 0.2 degrees (typical industrial-grade inclinometer accuracy) and the attitude distribution deviates significantly from the Gaussian assumption (i.e., the attitude parameter sequence has a p-value less than 0.05 after the Shapiro-Wilk normality test), the attitude error remaining from the first rotation in the preprocessing stage is transmitted to the inversion process through the observation data due to the frozen attitude parameters in the previous implementation method. This results in a velocity inversion deviation of approximately 6.0%, a defect location error of approximately 0.22 meters, and a defect type classification accuracy of approximately 86% after inversion convergence. This extended implementation method, by using the attitude parameter R(α) and model parameter m as inversion variables in joint inversion, reduces the attitude parameter estimation error from the initial 0.2 degrees to 0.04 degrees, the velocity inversion deviation to 2.5%, the defect location error to 0.13 meters, and the defect type classification accuracy to 94% upon inversion convergence. Regarding the inversion success rate, the previous implementation method achieved an inversion success rate of approximately 82% under these conditions, while this extended implementation method increases the success rate to over 93%. The aforementioned comparative experiments fully demonstrate that this preferred embodiment, through its paradigm innovation of attitude-waveform joint inversion, achieves a significant improvement in inversion accuracy under complex conditions where attitude measurement noise is high or the attitude distribution deviates from the Gaussian assumption. Furthermore, this preferred embodiment, along with the aforementioned three preferred embodiments based on neural networks, reinforcement learning, and self-supervised learning, as well as the extended embodiment of attitude-frequency coupled domain S-transform correction, can be used independently or in arbitrary combinations. When used in combination, the technical effects of each method can be superimposed to further improve inversion accuracy.
[0102] Step 4.5, Model parameter update and convergence criterion verification.
[0103] After constructing the adaptive update direction, the processor updates the model parameters according to the following expression: Where: m (k+1) This represents the updated model parameters for the (k+1)th iteration; m (k) Δm represents the model parameters in the current k-th iteration. (k) This indicates the amount of model updates in the current iteration.
[0104] The convergence criterion verification step is executed after each iteration when the model parameters are updated. The processor calculates the convergence coefficient ε for the current iteration. (k) The specific expression for the convergence coefficient is as follows: Where: ε (k) u represents the convergence coefficient of the current k-th iteration; obs This represents the observed wavefield elements at all time sampling points, consisting of all source numbers s and all detector numbers r. The observed wavefield tensor; Indicates in The forward wavefield tensor obtained by forward modeling the isotropic elastic wave equation has the following elements: When ε (k) When the convergence threshold is less than the preset convergence threshold (in this embodiment, the preset convergence threshold is 1×10), -4 Once the inversion iteration process is determined to be converged, the processor terminates the iteration process and uses the current model parameters as the final model parameters.
[0105] Step S5: Decoupling of transverse and longitudinal wave velocities and identification of pile defects.
[0106] After iterative convergence, the processor analyzes the final model parameters as three components: the spatial distribution λ of the first Lamé coefficient, the spatial distribution μ of the second Lamé coefficient, and the spatial distribution ρ of the density. Then, the processor calculates the three-dimensional distributions of the pile's shear wave velocity and longitudinal wave velocity, respectively. The pile's shear wave velocity v... s The specific expression is: Longitudinal wave velocity v of pile body p The specific expression is: Where: v s Indicates the shear wave velocity of the pile; v p λ represents the longitudinal wave velocity of the pile; λ and μ are both Lamé coefficients; ρ is the density.
[0107] The defect identification process, based on the pile design parameters (C35 concrete, design age 28 days), consults the concrete wave velocity reference range listed in the appendix of the "Technical Specification for Testing of Building Foundation Piles" (JGJ106-2014) to determine that the reference value for the longitudinal wave velocity of the intact concrete of this pile is approximately 4000 m / s and the reference value for the transverse wave velocity is approximately 2300 m / s. A preset deviation threshold of ±15% is set (i.e., the critical wave velocity reduction corresponding to a first-level reduction in concrete strength grade). The processor compares the inverted three-dimensional distribution of transverse wave velocity (vs) and the three-dimensional distribution of longitudinal wave velocity (vp) with the reference values at each spatial grid point. Spatial areas where the wave velocity deviation from the reference value exceeds the preset deviation threshold are identified as pile defect areas. Finally, a pile integrity status information report containing the spatial location (depth range, cross-sectional location) of the pile defect areas is output. (See [link to relevant documentation]). Figure 9 .
[0108] To demonstrate the unexpected technical advantages of the tilt-adaptive inversion proposed in this invention compared to the conventional full-waveform inversion without tilt modulation, a comparative experiment was conducted on the same pile foundation model in this embodiment. The results are as follows: Figure 6 As shown. Specific control experiment data are as follows: In the control group experiment without attitude confidence factor modulation (i.e., η(t) is always equal to 1), the full waveform inversion iteration still failed to reach the convergence threshold after 200 iterations, and the final objective function residual was still hovering around 26% of the initial value; in the experimental group with attitude confidence factor modulation, the full waveform inversion iteration reached the convergence threshold after only 80 iterations, and the objective function residual dropped to 0.8% of the initial value.
[0109] In terms of the accuracy of the inversion results, the control group had an error of 0.35 meters in the location of the pile diameter reduction defect (i.e., the difference between the identified defect depth and the actual defect depth) and failed to identify the diameter expansion defect; the experimental group had an error of only 0.12 meters in the location of the pile diameter reduction defect and successfully identified the diameter expansion defect with a location error of 0.18 meters.
[0110] In terms of noise robustness, when random noise with a signal-to-noise ratio of 10 dB is added to the observed signal, the control group inversion completely diverges, while the experimental group can still maintain stable convergence and keep the defect identification accuracy within 0.25 meters.
[0111] The above comparative experiments fully demonstrate that the present invention maps the physical quantity of attitude parameters to an attitude confidence factor η(t) using a Gaussian kernel and synchronously couples it to the tilt angle adaptive step size factor η. (k) The computational and gradient direction preconditioning stages achieve significant and unexpected technical improvements in convergence speed, inversion accuracy, and noise robustness compared to conventional full-waveform inversion without tilt modulation. Furthermore, by employing the aforementioned three preferred implementation methods based on neural networks and reinforcement learning, the stability and accuracy of the inversion can be further enhanced under more complex conditions (such as non-Gaussian attitude distributions, severe sensor noise, saddle point platforms, etc.).
[0112] The specific embodiments described in this specification are merely preferred embodiments of the present invention. Equivalent substitutions or obvious modifications made by those skilled in the art without departing from the spirit of the present invention—including but not limited to: using other types of borehole geophones instead of the three-component accelerometer, using other forms of coupling mechanisms instead of the coupling airbag, using other nonlinear mapping functions instead of the Gaussian kernel function, using other adaptive optimization algorithms instead of the specific calculation form of the tilt angle adaptive step size factor, and employing other network architectures to implement the aforementioned preferred embodiments of neural networks and reinforcement learning—should all be considered to fall within the protection scope of the present invention. σ The initial step size η0 is equal to 0.05, the adaptive adjustment coefficient γ is equal to 5.0, and the momentum term coefficient β is equal to 0.15 radians. (k) Equal to 0.9, convergence threshold 1×10 -4The values (such as ±15% deviation from the threshold) are all exemplary values, and those skilled in the art can make reasonable adjustments based on specific engineering geological conditions, pile foundation design parameters, detection accuracy requirements, and other factors.
[0113] The contents not described in detail in this specification are existing technologies known to those skilled in the art.
Claims
1. A method for detecting existing pile foundations based on borehole elastic waves, characterized in that, Includes the following steps: Elastic wave signals from the lateral seismic source of the pile foundation are sequentially collected at multiple observation nodes along the depth direction of the pile foundation by a probe unit lowered into the borehole, and the three-dimensional spatial attitude parameters and spatial position information of the probe are simultaneously acquired at each observation node. Based on the collected three-dimensional spatial attitude parameters, a coordinate system transformation is performed on the elastic wave signal to transform the elastic wave signal from the probe's body coordinate system to the absolute geographic reference coordinate system, thereby obtaining the processed elastic wave signal. The forward modeling wave field is obtained by performing forward modeling on an underground model containing pile foundation and surrounding medium based on the isotropic elastic wave equation. The underground model is characterized by model parameters. A full waveform inversion objective function is constructed based on the waveform difference between the forward-modeled simulated wavefield and the processed elastic wave signal. To minimize this objective function, a tilt-adaptive full waveform inversion iterative process is performed on the model parameters. Each iteration of this process corresponds to one of the multiple observation nodes, and each iteration includes the following steps: Attitude confidence factor calculation step: The tilt angle component, azimuth angle component and rotation angle component are extracted from the three-dimensional spatial attitude parameters of the observation node corresponding to the current iteration step to construct the tilt angle feature vector. With the pre-calibrated probe tilt-free ideal attitude as the center and the preset bandwidth parameter as the scale, the distance between the tilt angle feature vector and the tilt-free ideal attitude is nonlinearly mapped into a continuous scalar between zero and one as the attitude confidence factor using a Gaussian kernel function. The adaptive step size factor calculation step involves: taking the arithmetic mean of the attitude confidence factor over all time sampling points corresponding to the current observation node to obtain the attitude confidence average value for the current iteration; multiplying the preset adaptive adjustment coefficient by the square of the L2 norm of the difference between the tilt eigenvectors between two adjacent iterations, and using the negative of the product as the exponent of an exponential function with the natural constant as the base, to obtain the attitude change attenuation factor for the current iteration; and using the product of the preset initial step size, the attitude confidence average value, and the attitude change attenuation factor as the adaptive step size factor for the current iteration. Gradient direction calculation step: Based on the gradient vector of the full waveform inversion objective function with respect to the model parameters, the attitude confidence factor is used as a product factor to construct a precondition for the gradient vector to obtain the gradient direction of the current iteration; the adaptive step size factor is combined with the gradient direction to obtain the model update amount of the current iteration, and the model parameters are updated with the model update amount; Once the full waveform inversion iteration process satisfies the convergence criterion, the medium wave velocity distribution is calculated based on the model parameters at the time of convergence, and the pile integrity status or defect distribution is determined based on the medium wave velocity distribution.
2. The method according to claim 1, characterized in that, The adaptive adjustment coefficient is pre-calibrated based on the statistical characteristics of the difference in tilt feature vector between two adjacent iterations during the full waveform inversion iteration process, and is further configured as follows: when the attitude change amplitude represented by the statistical characteristics exceeds a preset threshold, the adaptive adjustment coefficient is increased to enhance the attenuation intensity of the attitude change attenuation factor; when the attitude change amplitude represented by the statistical characteristics is lower than the preset threshold, the adaptive adjustment coefficient is decreased to maintain the search capability.
3. The method according to claim 1, characterized in that, The gradient direction calculation step is specifically implemented through the following steps: calculating the gradient vector of the full waveform inversion objective function under the current model parameters based on the adjoint state method; The gradient vector is smoothed and corrected using a preprocessing operator, and the attitude confidence factor is used as a product factor in the smoothing and correction operation to obtain the corrected gradient vector. The momentum term is constructed by multiplying the preset momentum term coefficient by the model update amount of the previous iteration; The negative of the corrected gradient vector is added to the momentum term, and the result is used as the gradient direction for the current iteration.
4. The method according to claim 1, characterized in that, The process of performing coordinate system transformation on the elastic wave signal based on the acquired three-dimensional spatial attitude parameters includes: constructing a three-dimensional spatial rotation matrix based on the tilt, azimuth, and rotation components in the acquired three-dimensional spatial attitude parameters; performing an inverse spatial coordinate system transformation operation on the acquired elastic wave signal using the inverse matrix of the three-dimensional spatial rotation matrix, transforming the elastic wave signal from the probe body coordinate system to the absolute geographic reference coordinate system; and using the transformed elastic wave signal as the processed elastic wave signal.
5. The method according to claim 1, characterized in that, After each iteration of the full waveform inversion iterative process completes the update of the model parameters, the convergence criterion is checked to see if it is satisfied. The convergence criterion is that the convergence coefficient calculated based on the current iteration model parameters is less than a preset convergence threshold. The convergence coefficient is the ratio of the norm of the difference between the observed wavefield and the forward wavefield calculated based on the current iteration model parameters to the norm of the observed wavefield. The full waveform inversion iterative process is terminated when the convergence criterion is satisfied.
6. The method according to claim 1, characterized in that, The step of calculating the medium wave velocity distribution based on the model parameters at convergence includes: resolving the model parameters at convergence into three components: a first Lamé coefficient, a second Lamé coefficient, and density; calculating the shear wave velocity distribution of the pile body based on the arithmetic square root of the ratio of the second Lamé coefficient to the density; calculating the longitudinal wave velocity distribution of the pile body based on the arithmetic square root of the ratio of the sum of twice the first Lamé coefficient and the second Lamé coefficient to the density; and using the shear wave velocity distribution and the longitudinal wave velocity distribution as the medium wave velocity distribution.
7. The method according to claim 6, characterized in that, The step of determining the integrity status or defect distribution of the pile body based on the medium wave velocity distribution includes: determining a reference value for the wave velocity of intact concrete based on the pile body design parameters; comparing the shear wave velocity distribution of the pile body with the longitudinal wave velocity distribution of the pile body and the reference value for the wave velocity of intact concrete; identifying spatial areas where the wave velocity deviates from the reference value by more than a preset deviation threshold as pile body defect areas; and outputting pile body integrity status information including the spatial location of the pile body defect areas.
8. A device for detecting existing pile foundations based on borehole elastic waves, characterized in that, This includes probe units lowered into the borehole and data processing workstations located on the surface; The probe unit has a coupling airbag arranged in the circumferential direction on the outer wall. The probe unit integrates a three-component accelerometer for acquiring elastic wave signals, an inclinometer for measuring the three-dimensional spatial attitude parameters of the probe in real time, and a positioning device for recording the spatial position information of the probe. The three-component accelerometer, the tilt meter, and the positioning device share the same time reference and are configured to synchronously output the elastic wave signal, the three-dimensional spatial attitude parameters, and the spatial position information. The processor of the data processing workstation is connected to the probe unit via a real-time data link and is configured to execute the method as described in any one of claims 1 to 7; wherein the processor is further configured to: after the probe unit is lowered to each observation node, control the coupling airbag to inflate until the coupling airbag is in contact with the borehole wall; after the coupling airbag is in contact with the borehole wall, initiate the iteration step of the full waveform inversion iteration process corresponding to the observation node; within the iteration step, read the three-dimensional spatial attitude parameters output by the inclinometer at the observation node from the real-time data link, calculate the current attitude confidence factor based on the three-dimensional spatial attitude parameters, and send the attitude confidence factor to the adaptive step size factor calculation stage and the gradient direction calculation stage respectively; Before the processor completes the calculation of the model update amount for this iteration, it does not accept data from the next observation node of the real-time data link.
9. The apparatus according to claim 8, characterized in that, The real-time data link between the probe unit and the data processing workstation is implemented using optical fiber communication. The probe unit is further provided with a storage module, which is configured to temporarily store the elastic wave signal, the three-dimensional spatial attitude parameters and the spatial position information when the real-time data link is blocked, and transmit the temporarily stored data to the data processing workstation through the optical fiber communication after the real-time data link is restored.
10. An electronic device for pile foundation testing, characterized in that, The device includes at least one processor and a memory communicatively connected to the at least one processor, the memory storing instructions executable by the at least one processor; when executed by the at least one processor, the instructions cause the electronic device to perform the method as described in any one of claims 1 to 7.