Intelligent control method for TBM construction

Abstract

The application discloses a kind of TBM construction intelligent control method, belong to tunneling construction technical field.For the problem that control instruction is delayed due to the response delay of propulsion hydraulic system and the poor adaptability of fixed update strategy, the method collects cutterhead speed, torque, propulsion speed, propulsion force and body attitude angle deviation value and packs and stores as historical data set;Predict the propulsion force and torque adjustment amount by LSTM model;Valve core displacement and piston displacement signal are collected and stored in storage queue and arranged in reverse order as inversion response sequence;Determine the weight output valve core displacement compensation signal with average displacement deviation and adaptively update the queue;Valve core target displacement instruction and compensation signal are superimposed to drive electro-hydraulic proportional valve, and torque adjustment amount controls frequency converter.The method can compensate for response delay without increasing hardware, and improve adaptability.The method of the application is used for the intelligent control of TBM tunneling.

Description

Technical Field

[0001] This invention relates to the field of tunnel boring machine (TBM) construction technology, specifically to an intelligent control method for TBM construction. Background Technology

[0002] Full-face tunnel boring machines (TBMs) are increasingly used in the construction of long tunnels, but the setting of their tunneling control parameters still mainly relies on operators' subjective judgment and adjustments based on real-time monitoring data. Because the TBM tunneling process involves multi-parameter coupling and complex and variable surrounding rock conditions, it is difficult for operators with limited experience and energy to coordinate all information and provide control parameters that match the current rock mass conditions. This parameter selection method, which relies on manual experience, not only increases the workload but also frequently results in mismatches between control parameters and rock mass conditions during actual construction, hindering further improvements in tunneling efficiency.

[0003] To reduce reliance on human experience, existing technologies have proposed using data-driven models to predict tunneling parameters. Application CN112183993A discloses an intelligent prediction and optimization decision-making method for TBM tunneling control parameters. This method preprocesses TBM tunneling parameters and cutterhead vibration acceleration, uses a long short-time memory network model to predict the tunneling control parameters, and generates optimized tunneling control parameters based on the prediction results. Application CN110895730A discloses a TBM tunneling parameter prediction method based on the LSTM algorithm. This method uses cutterhead rotation speed, feed rate, total thrust, and cutterhead torque as training features to establish a prediction model, outputting predicted values ​​of tunneling parameters in the stable phase, providing operators with suggestions for subsequent operating parameters.

[0004] However, the above-mentioned scheme uses the predicted tunneling parameters for subsequent control or operational decisions without compensating for the inherent response delay between the control command received by the valve core of the electro-hydraulic proportional valve and the corresponding displacement of the piston in the propulsion hydraulic cylinder. When the surrounding rock conditions change, the actual applied propulsion force lags behind the change in control command, and the control command cannot keep up with the changes in working conditions. In addition, existing technologies mostly use fixed-cycle or single-criterion triggering for online model updates, which makes it difficult to balance data stability with the ability to track dynamic changes in tunneling status.

[0005] The TBM propulsion hydraulic system is the core actuator driving the cutterhead forward to break rock, consisting of multiple hydraulic cylinders operating in parallel. During tunneling, there is an inherent response delay between the receipt of a control command by the valve core of the electro-hydraulic proportional valve and the corresponding displacement response of the piston in the propulsion hydraulic cylinder. This delay is determined by physical factors such as the compressibility of the hydraulic oil, the friction resistance along the pipeline, and the hydraulic resistance at the valve orifice. Slow response during the initial startup phase reduces the system's initial tracking performance and limits overall control accuracy. Under complex conditions such as alternating layers of soft and hard rock and high water content, this delay can easily induce problems such as speed lag, overshoot, or pressure spikes, thereby amplifying synchronization errors and leading to localized load concentration. During steady-state tunneling, when the surrounding rock conditions change, the response delay of the propulsion hydraulic system causes the actual applied propulsion force to lag behind the changes in control commands, resulting in control commands failing to keep pace with changes in the working conditions.

[0006] The response characteristics of electro-hydraulic proportional valves are affected by physical parameters such as hydraulic oil temperature and viscosity, while the dynamic characteristics of the propulsion hydraulic cylinder depend on load conditions and system pressure. These factors mean that the specific duration of the response delay is not constant during tunneling. Enabling the actual output of the propulsion hydraulic system to follow changes in control commands more promptly, without requiring additional hardware investment, is an engineering challenge in the field of TBM automatic control. Summary of the Invention

[0007] This invention provides an intelligent control method for TBM construction, which compensates for the response delay of the propulsion hydraulic system without increasing hardware investment, so that the actual output of the propulsion hydraulic system can follow the changes in control commands more promptly, thereby improving the adaptability of the control process to changes in surrounding rock conditions.

[0008] To achieve these objectives and other advantages according to the present invention, a smart control method for TBM construction is provided, comprising the following steps:

[0009] Step 1: Collect cutterhead rotation speed, cutterhead torque, propulsion speed, propulsion force, and machine attitude angle deviation values ​​at a fixed sampling period, and package them into data frames and store them in the historical tunneling dataset;

[0010] Step 2: Input the propulsion speed, mean propulsion force, variance of propulsion force, fuselage attitude angle deviation and rate of change of the historical tunneling dataset from 10 to 50 data frames before the current time step into the long short-term memory network model to obtain the propulsion force adjustment and torque adjustment.

[0011] Step 3: Collect valve core displacement signal and piston displacement signal and store them in historical tunneling dataset and storage queue. The length of the dataset is equal to the number of sampling cycles from receiving the command from the valve core to the piston producing displacement. After the dataset is full, arrange them in reverse order to form the inversion response sequence.

[0012] Step 4: Using a preset number of rotations of the cutterhead as the judgment period, calculate the average displacement deviation between the actual piston displacement value within this judgment period and the historical piston displacement value in the historical tunneling data one delay time earlier. The ratio of this average displacement deviation to the preset maximum deviation is used as the second weight. When the value of the second weight is greater than 1, it is taken as 1. 1 minus the second weight is used as the first weight. The current valve core displacement signal and the valve core displacement signal at the corresponding time in the inversion response sequence are weighted and summed according to the two weights to obtain the valve core displacement compensation signal. When the average displacement deviation is less than the first deviation threshold, the storage queue is updated with the first update ratio; when it is greater than the second deviation threshold, the storage queue is updated with the second update ratio. Other cases remain unchanged. The second deviation threshold is greater than the first deviation threshold.

[0013] Step 5: Convert the thrust adjustment amount into a valve core target displacement command and superimpose it with the compensation signal before sending it to the electro-hydraulic proportional valve; convert the torque adjustment amount into a torque control signal and send it to the frequency converter of the cutter head drive system.

[0014] Preferably, in the TBM construction intelligent control method, in step one, the current signal of the cutterhead drive motor and the pressure signal of the propulsion hydraulic system are also collected simultaneously; between step one and step two, the following steps are also included, wherein one collection cycle is defined as one revolution of the cutterhead:

[0015] Envelope detection is performed on the current signal acquired in step one to obtain the current amplitude envelope signal. Fourier transform is performed on the current amplitude envelope signal to obtain the envelope spectrum. The fundamental frequency amplitude and the second harmonic amplitude of the cutter head rotation frequency are extracted from the envelope spectrum. The ratio of the second harmonic amplitude to the fundamental frequency amplitude is calculated and recorded as the second harmonic ratio.

[0016] For the pressure signal of the propulsion hydraulic system collected in step one, calculate the number of rising edges of the pressure signal per unit time, and record the number of rising edges as the pressure impact frequency.

[0017] The second harmonic ratio is multiplied by the pressure impact frequency to obtain the geological feature fusion value;

[0018] When the geological feature fusion value increases monotonically over three consecutive acquisition cycles, it is determined that the rock at the working face has hardened. The mean thrust extracted in step two is multiplied by a first correction coefficient greater than 1 and then used as the corrected mean thrust input into the long short-term memory network model.

[0019] When the geological feature fusion value decreases monotonically over three consecutive acquisition cycles, the lithology of the tunnel face is determined to be softened. The mean thrust extracted in step two is multiplied by a second correction coefficient less than 1 and then used as the corrected mean thrust input into the long short-term memory network model.

[0020] Preferably, in the TBM construction intelligent control method, the current signal of the cutterhead drive motor and the pressure signal of the propulsion hydraulic system are collected. The current signal collected in step one is subjected to envelope detection to obtain the current amplitude envelope signal. The specific method for obtaining the envelope spectrum by performing Fourier transform on the current amplitude envelope signal is as follows:

[0021] One acquisition cycle is defined as one revolution of the cutter head, and current and pressure signals are continuously acquired within each acquisition cycle.

[0022] The current signal of the current acquisition cycle and the previous four consecutive acquisition cycles, a total of five acquisition cycles, are stored in the current sliding window; the pressure signal of the current acquisition cycle and the previous two consecutive acquisition cycles, a total of three acquisition cycles, are stored in the pressure sliding window.

[0023] At the end of each acquisition cycle, the envelope spectrum is obtained by performing envelope detection and Fourier transform on all current signals within the current sliding window. The number of rising edges is counted for all pressure signals within the pressure sliding window. This number is divided by the total duration of the pressure sliding window (in minutes) to obtain the pressure impact frequency of the pressure sliding window.

[0024] The geological feature fusion value for that acquisition period is obtained by multiplying the second harmonic ratio of the current sliding window with the pressure impact frequency of the pressure sliding window.

[0025] Preferably, in the TBM construction intelligent control method, at the end of each acquisition cycle, the envelope spectrum is obtained by performing envelope detection and Fourier transform on all current signals within the current sliding window, and the number of rising edges is counted for all pressure signals within the pressure sliding window. This number is then divided by the total duration of the pressure sliding window (in minutes) to obtain the pressure impact frequency of the pressure sliding window.

[0026] At the end of each acquisition cycle, the intersection of the current sliding window and the pressure sliding window on the time axis is defined as the common time window, and the length of the common time window is equal to the length of the pressure sliding window.

[0027] After envelope detection of all current signals within the current sliding window, Fourier transform is performed to obtain the envelope spectrum, and the second harmonic ratio is extracted from the envelope spectrum.

[0028] The number of rising edges of the pressure signal within the common time window is counted, and this number is divided by the total duration of the common time window (in minutes) to obtain the pressure impact frequency of the common time window.

[0029] The geological feature fusion value for that acquisition period is obtained by multiplying the second harmonic ratio by the pressure impact frequency of the common time window.

[0030] Preferably, in the TBM construction intelligent control method, the specific method for obtaining the thrust adjustment and torque adjustment from the long short-term memory network model in step two is as follows:

[0031] The Long Short-Term Memory network model outputs a thrust prediction and a torque prediction at each time step;

[0032] The instantaneous thrust error is obtained by subtracting the predicted thrust value output at the current time step from the actual thrust value collected at the current time step, and the instantaneous torque error is obtained by subtracting the predicted torque value output at the current time step from the actual torque value collected at the current time step.

[0033] The instantaneous thrust errors of P consecutive time steps, including the current time step, are used to form a thrust error sequence. The variance of the thrust error sequence is calculated and denoted as the thrust error variance. The instantaneous torque errors of P consecutive time steps, including the current time step, are used to form a torque error sequence. The variance of the torque error sequence is calculated and denoted as the torque error variance.

[0034] When the thrust error variance is greater than the first variance threshold, the thrust prediction value of the current time step is replaced with the arithmetic mean or weighted moving average of the actual thrust values ​​of the previous preset number of time steps and then output as the thrust adjustment amount. The preset number is less than P.

[0035] When the torque error variance is greater than the second variance threshold, the torque prediction value of the current time step is replaced with the arithmetic mean or weighted moving average of the actual torque values ​​of the previous preset number of time steps, and then used as the torque adjustment amount output.

[0036] Preferably, in the TBM construction intelligent control method, when the propulsion error variance is greater than the first variance threshold, the specific method for replacing the predicted propulsion value of the current time step with the arithmetic mean or weighted moving average of the actual propulsion values ​​of a preset number of time steps prior to the current time step as the output of the propulsion adjustment amount is as follows:

[0037] Calculate the variance of the thrust error sequence for P consecutive time steps, including the current time step, which is the thrust error variance, denoted as the current variance;

[0038] For each of the P consecutive time steps including the current time step, calculate the difference between the predicted thrust value of the time step and the actual thrust value offset by k time steps before the time step, where k ranges from 1 to K and K is the preset maximum delay step number, and obtain K offset error sequences. Calculate the variance of each offset error sequence and denote it as the offset variance.

[0039] The current variance is compared with each offset variance one by one. When there is an offset variance smaller than the current variance, the k value corresponding to the smallest offset variance is selected as the estimated delay step; when each offset variance is greater than or equal to the current variance, the estimated delay step is 0.

[0040] When the estimated delay step is 0 or 1, the predicted thrust value of the current time step is replaced with the actual thrust value of the previous time step and then output as the thrust adjustment amount.

[0041] When the estimated delay step number is greater than 1 and less than or equal to the preset first delay threshold, the predicted thrust value of the current time step is replaced with the arithmetic mean of the actual thrust values ​​of the estimated delay step number of time steps before the current time step, and then the thrust adjustment amount is output.

[0042] When the estimated delay step count is greater than the first delay threshold, the predicted thrust value at the current time step is replaced with the weighted moving average of the actual thrust values ​​at a preset fixed number of time steps prior to the current time step, and this is used as the thrust adjustment amount output.

[0043] Preferably, in the TBM construction intelligent control method, the specific method for calculating the variance of the thrust error sequence over P consecutive time steps, including the current time step, and the specific method for calculating the variance of each offset error sequence is as follows:

[0044] Perform linear regression on the thrust error sequence to obtain the linear regression slope of the thrust error sequence;

[0045] A trend term sequence is generated based on the slope and intercept of the linear regression. The trend value at the corresponding time in the trend term sequence is subtracted from each value in the thrust error sequence to obtain the detrended thrust error sequence.

[0046] Calculate the variance of the detrended driving force error sequence and use it as the current variance;

[0047] After performing linear detrending on each offset error sequence, the variance of each offset error sequence is calculated and used as the offset variance.

[0048] Preferably, in the TBM construction intelligent control method, the long short-term memory network model in step two is updated online in the following manner:

[0049] An update cycle is defined as N revolutions of the cutterhead, where N is an integer greater than 1. At the end of each update cycle, the propulsion prediction sequence and torque prediction sequence output by the long short-term memory network model at each time step within that update cycle are collected, along with the actual propulsion sequence and actual torque sequence at the corresponding time.

[0050] Calculate the difference between the predicted values ​​of two adjacent time steps in the thrust prediction value sequence to obtain the thrust prediction difference sequence. Calculate the difference between the actual values ​​of two adjacent time steps in the actual thrust value sequence to obtain the actual thrust difference sequence. Calculate the sign consistency rate between the thrust prediction difference sequence and the actual thrust difference sequence, denoted as the thrust trend consistency rate.

[0051] Calculate the difference between the predicted values ​​of two adjacent time steps in the torque prediction value sequence to obtain the torque prediction difference sequence. Calculate the difference between the actual values ​​of two adjacent time steps in the actual torque value sequence to obtain the actual torque difference sequence. Calculate the sign consistency rate between the torque prediction difference sequence and the actual torque difference sequence, denoted as the torque trend consistency rate.

[0052] When the thrust trend consistency rate is lower than the first consistency rate threshold or the torque trend consistency rate is lower than the second consistency rate threshold, all data frames within the update period are added as new training samples to the training set of the long short-term memory network model. The model parameters are then updated after incremental training of the long short-term memory network model.

[0053] Preferably, the TBM construction intelligent control method further includes the following steps:

[0054] When the thrust trend consistency rate is greater than or equal to the first consistency rate threshold, the relative deviation between the thrust prediction value and the actual thrust value at each time step in the update cycle is further calculated. The relative deviation is the thrust prediction value minus the actual thrust value and divided by the non-zero actual thrust value to obtain the thrust relative deviation sequence.

[0055] Calculate the slope of the linear regression of the relative deviation sequence of propulsion force;

[0056] When the slope of the linear regression is positive and greater than the preset positive slope threshold, or when the slope of the linear regression is negative and the absolute value is greater than the absolute value of the preset negative slope threshold, all data frames in the update period are added as new training samples to the training set of the long short-term memory network model, and the model parameters are updated after incremental training of the long short-term memory network model.

[0057] Preferably, the TBM construction intelligent control method further detects abrupt amplitude drift in the following manner:

[0058] Calculate the relative deviation between the predicted thrust value and the actual thrust value at each time step within the update cycle, arrange them in chronological order to obtain the thrust relative deviation sequence, and divide the sequence into a first segment sequence and a second segment sequence in chronological order. When the sequence length is odd, the middle element is included in the first segment sequence.

[0059] Calculate the arithmetic mean of the first part of the sequence as the first part mean, and calculate the arithmetic mean of the second part of the sequence as the second part mean.

[0060] Calculate the difference between the mean of the second segment and the mean of the first segment;

[0061] When the difference is greater than the first difference threshold, or when the difference is less than the second difference threshold, all data frames in the update period are added as new training samples to the training set of the long short-term memory network model. The model parameters are updated after incremental training of the long short-term memory network model. The second difference threshold is negative.

[0062] The present invention has at least the following beneficial effects:

[0063] 1. This invention uses a weighted superposition of the valve core displacement signal at the current time step and the valve core displacement signal at the corresponding moment in the inverted response sequence as the output compensation signal. This integrates historical response information into the drive signal of the electro-hydraulic proportional valve, improving the response characteristics of the propulsion hydraulic system without increasing hardware costs. The update ratio of the storage queue is adaptively adjusted by the magnitude of the average displacement deviation. When the average displacement deviation is small, a low update ratio is used to keep the queue stable; when the average displacement deviation is large, a high update ratio is used to enable the queue to quickly track system dynamic changes, thus balancing queue stability and adaptability.

[0064] 2. This invention calculates the second harmonic ratio by collecting the current signal of the cutterhead drive motor. The ratio of the second harmonic to the fundamental frequency amplitude in the envelope spectrum reflects the degree of nonlinear distortion of the rock-breaking load waveform, thus avoiding the need for additional vibration sensors to obtain geological information. By collecting the pressure signal of the propulsion hydraulic system and statistically analyzing the pressure impact frequency, the impact frequency of the rock-breaking load is reflected in the time domain. Multiplying the frequency domain features with the time domain features yields a fused geological feature value. The monotonic trend of this fused value over multiple consecutive acquisition cycles indicates the gradual change direction of lithology. Based on this, the input values ​​of the long short-term memory network model are corrected, helping to correct the model input before it deviates from the actual geological conditions, thereby reducing the occurrence of mispredictions.

[0065] 3. This invention sets a sliding window of five acquisition cycles for the current signal, which helps ensure that the signal length participating in the Fourier transform is sufficient to obtain the frequency resolution required to distinguish the fundamental frequency and the second harmonic. Setting a sliding window of three acquisition cycles for the pressure signal avoids over-smoothing of transient impact information due to excessively long windows. The end times of both windows are aligned with the end of the current acquisition cycle, ensuring that the second harmonic ratio and the pressure impact frequency are naturally synchronized in time. The lengths of the two windows are respectively matched to the time scale requirements of their respective feature extraction methods, which helps improve the accuracy of feature values.

[0066] 4. This invention defines the intersection of the current sliding window and the pressure sliding window on the time axis as a common time window, so that the statistics of pressure impact frequency are strictly based on the time period covered by both windows. This can eliminate the interference of the two additional acquisition cycles in the current sliding window compared to the pressure sliding window on time alignment, reduce the problem of incomplete overlap of time periods caused by different starting times of sliding windows of different lengths, and make the two multiplied feature values ​​have a more accurate correspondence on the time axis, which is beneficial to improving the accuracy of geological feature fusion values.

[0067] 5. This invention obtains the instantaneous error by subtracting the predicted value from the actual value at the current time step. It then constructs an error sequence from the instantaneous errors of the current time step and several consecutive previous time steps and calculates the variance. The variance, reflecting the overall fluctuation of the error sequence, is used to determine the reliability of the current predicted value. When the error variance exceeds a threshold, the current predicted value is deemed unreliable and discarded. Instead, the arithmetic mean of the actual values ​​within a shorter time window preceding the current time step is used as the adjustment output. This ensures that control commands are generated based on the latest actual tunneling status even when the model generates delayed predictions, helping to reduce the adverse effects of prediction delays on control performance. Furthermore, this method does not require altering the structure of the Long Short-Term Memory network model itself, has a low computational load, and is easily implemented on existing industrial control platforms.

[0068] 6. This invention compares the current variance with each offset variance one by one, and estimates the current prediction delay steps using the principle of minimizing the variance after offset. A tiered processing strategy is adopted based on the estimated delay steps: when the delay is no more than one sampling period, the actual value of the latest time step is used to replace the predicted value; when the delay step is small, the arithmetic mean of the actual values ​​of the time steps corresponding to the estimated delay step is used to replace the predicted value; when the delay step is large, the weighted moving average of the actual values ​​of a preset fixed number of time steps is used to replace the predicted value. This ensures that the compensation method matches the degree of delay, helping to avoid over-intervention in mild delays and under-intervention in severe delays.

[0069] 7. This invention performs linear regression on the propulsion error sequence to obtain the slope and intercept. Based on the slope and intercept, a trend term sequence is generated, and the trend term is subtracted from the error sequence. This eliminates the interference of trend components caused by factors such as tool wear or slow changes in the surrounding rock on the calculation of the error sequence variance. The detrended variance better reflects the degree of fluctuation of the error sequence around its mean, making the comparison results between the current variance and each offset variance more accurately reflect the delay step count, which helps improve the reliability of the delay step count estimation.

[0070] 8. In this invention, the sign consistency rate between the predicted difference sequence and the actual difference sequence is calculated. The trend consistency rate is used to quantitatively evaluate whether the direction of change of the predicted value is consistent with the actual value. When the trend consistency rate is lower than a threshold, it indicates that the direction predicted by the model has begun to deviate from the actual direction, at which point incremental training is triggered. The trend consistency rate is not sensitive to the dimensions and magnitude, and is less affected by changes in the magnitude of propulsion force and torque, which helps to trigger updates in the early stages when the model's predicted trend begins to deviate from the actual trend. The update cycle is based on the number of rotations of the cutterhead rather than a fixed time, so that the model update rhythm is naturally synchronized with the changing rhythm of the tunneling conditions.

[0071] 9. This invention, under the premise of achieving a consistent trend rate, further calculates the linear regression slope of the relative deviation sequence. The sign and magnitude of the slope are used to determine whether there is a persistent amplitude drift between the predicted and actual values. When the slope is positive and greater than a positive slope threshold, a positive amplitude drift is identified; when the slope is negative and its absolute value is greater than a negative slope threshold, a negative amplitude drift is identified. This allows the model update to simultaneously respond to both directional deviation and amplitude drift degradation modes. Even when the directional prediction is accurate but the amplitude continues to deviate, model updates can still be triggered, which helps avoid situations where the model is only updated after the amplitude drift has accumulated to a significant extent.

[0072] 10. This invention divides the relative deviation sequence into a first segment and a second segment in chronological order. The arithmetic mean of each segment is calculated, and the difference is then taken. The difference between the first and second segment means is used to detect abrupt amplitude drift. This method, along with the linear regression slope method, detects amplitude drift from different perspectives, thus complementing each other. When amplitude drift occurs near the middle of the update cycle, the difference between the two segment means can be reflected quickly, enabling the model update mechanism to respond simultaneously to both gradual and abrupt amplitude drift modes. It can trigger model updates earlier when abrupt amplitude drift occurs, which helps to shorten the response time from the occurrence of amplitude drift to model update.

[0073] Other advantages, objectives and features of the present invention will become apparent in part from the following description, and in part from those skilled in the art through study and practice of the invention. Detailed Implementation

[0074] The present invention will now be described in further detail so that those skilled in the art can implement it based on the description.

[0075] It should be understood that terms such as “having,” “comprising,” and “including” as used herein do not exclude the presence or addition of one or more other elements or combinations thereof.

[0076] According to one embodiment of the present invention, a TBM construction intelligent control method is provided. During the tunneling process, five parameters are collected in real time by a data acquisition module installed in the TBM main control box: cutterhead rotation speed, cutterhead torque, propulsion speed, propulsion force, and body attitude angle deviation. The body attitude angle deviation is determined by the deviation angle of the body attitude angle relative to the design axis, collected by an inclination sensor installed on the TBM body, and is expressed in degrees. The unit of cutterhead rotation speed is revolutions per minute (rpm), the unit of cutterhead torque is kilonewton-meters (kNm), the unit of propulsion speed is millimeters per minute (mm / min), and the unit of propulsion force is kilonewtons (kN). The sampling period of the data acquisition module can be set to 1 second, or any value between 0.5 seconds and 2 seconds as needed. In this invention, the sampling period is the smallest time unit for data acquisition and processing, with each time step corresponding to one sampling period. The acquisition period is calculated per revolution of the cutterhead, while the judgment period and update period are calculated per preset number of revolutions and N revolutions of the cutterhead, respectively. The duration of each of the acquisition period, judgment period, and update period consists of several sampling periods. At the end of each sampling cycle, the five data points collected in the current cycle—cutterhead rotation speed, cutterhead torque, propulsion speed, propulsion force, and fuselage attitude angle deviation—are packaged into a single data frame according to a fixed data frame format. This data frame is then stored sequentially in the historical tunneling dataset. The historical tunneling dataset is stored on the hard drive or solid-state storage of the TBM industrial control computer for data retrieval in subsequent steps.

[0077] The Long Short-Term Memory (LSTM) network model is deployed on the processor of the TBM industrial control computer. During prediction, 10 to 50 data frames prior to the current time window are read from the historical tunneling dataset. The specific window length can be preset within this range according to the actual tunneling conditions; for example, a window length of 30 data frames can be set in homogeneous strata, while a window length of 15 data frames can be set in strata with frequent lithological changes. Five feature values ​​are extracted from these data frames: thrust velocity, mean thrust force, thrust force variance, fuselage attitude angle deviation, and fuselage attitude angle deviation rate of change. The mean thrust force is the arithmetic mean of the thrust forces of each data frame within the window; the thrust force variance is the average of the squares of the differences between the thrust forces of each data frame within the window and the mean thrust force; the fuselage attitude angle deviation is the deviation angle of the current fuselage attitude angle relative to the design axis; and the fuselage attitude angle deviation rate of change is the difference between the current attitude angle deviation and the attitude angle deviation value at the previous time step. These five feature values ​​are organized into an input vector according to the format requirements of the LSM network model's input layer and input into the pre-trained LSM network model. After forward computation, the Long Short-Term Memory (LSTM) network model outputs predicted thrust and torque values ​​from its output layer, from which thrust and torque adjustments are derived. The unit for thrust adjustment is kilonewtons (kN), and the unit for torque adjustment is kilonewton-meters (kNm).

[0078] A valve core displacement sensor is installed on the electro-hydraulic proportional valve of the propulsion hydraulic system. This sensor can be a linear variable differential transformer type displacement sensor or a magnetostrictive displacement sensor, used to measure the axial displacement of the valve core in real time. A piston displacement sensor is installed on the piston rod of the propulsion hydraulic cylinder. This sensor can be a wire type displacement sensor or a magnetic grating type displacement sensor, used to measure the extension length of the piston rod in real time. The valve core displacement signal and the piston displacement signal are synchronously acquired with the same sampling period as in step one, and stored sequentially in a storage queue according to time order. The storage queue is a first-in, first-out buffer in the TBM industrial control computer memory, and its length is equal to the total number of sampling periods experienced by the propulsion hydraulic system from the moment the valve core of the electro-hydraulic proportional valve receives the control command until the piston of the propulsion hydraulic cylinder produces a displacement response corresponding to the command. The product of the total number of sampling periods and the sampling period is the delay time of the propulsion hydraulic system. The total number of sampling periods can be predetermined through a step response test: a step control signal is sent to the electro-hydraulic proportional valve, and the response curves of the valve core displacement and piston displacement are acquired and recorded at a frequency of not less than 100 Hz. The time difference between the moment the control signal is issued and the moment when the piston displacement reaches a specified percentage of its steady-state value, such as 90% of the steady-state value, is divided by the sampling period to obtain the length of the storage queue. When the hydraulic oil temperature or load conditions change significantly, the above step response test can be repeated to correct the length of the storage queue. When the amount of data stored in the storage queue reaches its length, all the data in the queue are arranged in reverse order according to the acquisition time from latest to earliest, forming an inversion response sequence. Simultaneously with the storage queue, the piston displacement signal is also stored in the historical tunneling dataset in chronological order, used to determine the actual value of the piston displacement within the cycle and the reading of its historical values.

[0079] A judgment cycle is triggered after each preset number of cutterhead rotations. The preset number of rotations can be set from 1 to 5 rotations; for example, it can be set to 3 rotations in homogeneous strata and 1 rotation in strata with frequent lithological changes. At each time step within the judgment cycle, the actual value of the hydraulic cylinder piston displacement at that time step is read from the historical tunneling dataset. Simultaneously, the historical piston displacement value at a historical moment delayed by one time interval is read from the historical tunneling dataset. The difference between the two values ​​is subtracted, and the absolute value is taken to obtain the displacement deviation value for that time step. At the end of the judgment cycle, the displacement deviation values ​​for all time steps within the judgment cycle are summed and divided by the total number of time steps in the judgment cycle to obtain the average displacement deviation. The unit of the average displacement deviation is millimeters.

[0080] The preset maximum deviation is the maximum permissible deviation of piston displacement pre-set according to the design parameters of the propulsion hydraulic system. Its value can be set within the range of 0.5 mm to 5 mm, depending on the accuracy requirements of the specific hydraulic system; for example, it can be set to 2 mm. The average displacement deviation is divided by the preset maximum deviation to obtain the ratio, which is used as the second weight; when the value of the second weight is greater than 1, it is set to 1. The first weight is equal to 1 minus the second weight. The larger the average displacement deviation, the more significant the deviation between the actual piston displacement value and the historical value before the delay. In this case, the second weight increases, that is, the proportion of the historical displacement signal increases, which is beneficial for strengthening the anticipatory compensation for the response delay. The smaller the average displacement deviation, the closer the actual displacement is to the historical displacement. In this case, the first weight increases, that is, the proportion of the current displacement signal increases, which is beneficial for maintaining the accuracy of control. At each time step, the valve core displacement signal at the current moment is read from the storage queue, and the valve core displacement signal at the historical moment corresponding to the current time step is read from the inversion response sequence. The valve core displacement signal at the current moment is multiplied by the first weight, and the valve core displacement signal at the historical moment is multiplied by the second weight. The two products are added together, and the result is output as the valve core displacement compensation signal.

[0081] The update ratio of the storage queue is determined based on the magnitude of the average displacement deviation. The first and second deviation thresholds are determined based on historical displacement deviation statistics under normal operating conditions of the propulsion hydraulic system. The first deviation threshold can be set within the range of 0.5 to 1.0 times the average historical displacement deviation, for example, it can be set to 0.8 times the average historical displacement deviation. The second deviation threshold can be set within the range of 1.5 to 3.0 times the average historical displacement deviation, for example, it can be set to 2.0 times the average historical displacement deviation. When the average displacement deviation is less than the first deviation threshold, the current tunneling state is determined to be stable. The valve core displacement signal and piston displacement signal collected within the current judgment period are updated in the storage queue according to the first update ratio. That is, the data stored in the storage queue is replaced according to the first update ratio. The amount of data replaced is the product of the storage queue length and the update ratio. If it is less than one sampling period, it is counted as one sampling period. The earliest data of equal quantity replaced is removed. The first update ratio can be set within the range of 10% to 30%, for example, it can be set to 20%. When the average displacement deviation exceeds the second deviation threshold, the current tunneling state is determined to be unstable. The valve core displacement signal and piston displacement signal collected during the current judgment period are then updated to the storage queue in the same manner using the second update ratio. The second update ratio can be set within the range of 70% to 100%, for example, it can be set to 80%. The second deviation threshold is greater than the first deviation threshold. When the average displacement deviation is greater than or equal to the first deviation threshold and less than or equal to the second deviation threshold, the current tunneling state is determined to be in a transition range, and the contents of the storage queue remain unchanged.

[0082] The thrust adjustment output from the Long Short-Term Memory (LSTM) network model is converted into an analog voltage or current signal by a digital-to-analog converter (DAC) module and then output as a thrust control signal. This thrust control signal is the target displacement command of the electro-hydraulic proportional valve spool, measured in millimeters. In this invention, the thrust adjustment and torque adjustment refer to the thrust and torque values ​​used to generate the control command after the following corrections, rather than increments relative to the previous moment. The target displacement command of the spool is added to the spool displacement compensation signal output in step four, and the result is sent as the final drive signal to the power amplifier of the electro-hydraulic proportional valve in the propulsion hydraulic system. This drives the spool to move to the target position, thereby adjusting the flow rate of oil entering the propulsion hydraulic cylinder and thus adjusting the thrust. The torque adjustment output from the LSM model is converted into an analog voltage or current signal by a DAC module and then sent as a torque control signal to the frequency converter of the cutterhead drive system. The frequency converter adjusts its output frequency and output voltage according to the received torque control signal, thereby adjusting the output torque of the cutterhead drive motor and thus adjusting the cutterhead torque.

[0083] The above steps run continuously in a loop during execution: after step five completes the control signal transmission, it returns to step one to continue collecting data for the next sampling period, forming a continuous control closed loop.

[0084] Compared to existing technologies that rely on operators' subjective experience to set tunneling parameters or use open-loop control methods that are difficult to cope with changes in surrounding rock conditions, this invention uses a long short-term memory network model to predict tunneling parameters. At the same time, it uses a storage queue and inversion response sequence to compensate the drive signal of the electro-hydraulic proportional valve on a sampling cycle basis, and adaptively adjusts the update ratio of the storage queue through the average displacement deviation. This helps to make the response of the propulsion hydraulic system more timely and accurate, and improves the adaptability of the control process to changes in surrounding rock conditions.

[0085] According to another embodiment of the present invention, in the provided intelligent control method for TBM construction, based on the synchronous acquisition of current and pressure signals in step one, a geological sensing and input correction stage is added between step one and step two. The cutterhead drive motor is generally an AC asynchronous motor or a permanent magnet synchronous motor, with a current transformer or Hall effect current sensor installed on its stator windings to acquire the three-phase current signal of the motor. The current sensor can be a closed-loop Hall effect current sensor, whose output signal is an analog voltage or current signal proportional to the primary current. The range is selected according to the rated current of the motor, generally 1.2 to 2.0 times the rated current. After being filtered by the signal conditioning circuit, the current signal is sent to the data acquisition card built into the industrial control computer. The sampling frequency can be set to 1 kHz to 10 kHz to ensure that the cutterhead rotation frequency and its integer multiples in the current amplitude envelope signal are completely preserved. The pressure signal of the propulsion hydraulic system is acquired by a pressure sensor installed near the oil inlet of the propulsion hydraulic cylinder. The pressure sensor can be a piezoresistive or strain gauge pressure sensor, and its range is selected according to the rated working pressure of the propulsion hydraulic system. Generally, it is selected as 1.5 to 2.0 times the rated pressure. The output signal is a current signal of 4 mA to 20 mA or a voltage signal of 0 V to 10 V. After signal conditioning, it is sent to the data acquisition card.

[0086] For the acquired current signal, envelope detection is first performed to obtain the current amplitude envelope signal. Envelope detection can be achieved using Hilbert transform or rectifier filtering. The carrier frequency of the stator current is the motor power supply frequency. The load fluctuation caused by the cutterhead rotation is reflected in the current amplitude envelope signal through amplitude modulation. Therefore, spectral analysis of the envelope signal can extract the cutterhead rotation frequency and its second harmonic component. Then, the current amplitude envelope signal is de-DC biased, and a Hanning window or Hamming window is applied to reduce spectral leakage. Finally, a fast Fourier transform is performed to obtain the envelope spectrum. The frequency resolution of the envelope spectrum depends on the duration of the signal participating in the Fourier transform; the longer the signal duration, the higher the frequency resolution. In the envelope spectrum, the fundamental frequency of the cutterhead rotation is the cutterhead rotation speed divided by 60, in Hertz, and its second harmonic is twice the fundamental frequency. The above feature extraction is performed under tunneling conditions, i.e., when the cutterhead rotation speed is greater than zero. The amplitude at the fundamental frequency and the amplitude at the second harmonic frequency are extracted from the envelope spectrum, both in amperes. The ratio of the second harmonic amplitude to the fundamental frequency amplitude is calculated. This dimensionless ratio is denoted as the second harmonic ratio. This ratio reflects the nonlinearity of the cutterhead's rock-breaking load: when the cutterhead cuts hard rock, the load fluctuates violently, the current waveform undergoes nonlinear distortion, and the proportion of the second harmonic component relative to the fundamental frequency component increases; when cutting soft rock, the load is stable, the current waveform is close to sinusoidal, and this ratio is smaller. Since this ratio is the ratio of the amplitudes of two frequency components in the same envelope spectrum, the influence of the overall amplitude change of the current signal on this ratio is eliminated, making this characteristic quantity mainly reflect the degree of distortion of the load waveform.

[0087] The acquired pressure signal is first subjected to low-pass filtering to eliminate high-frequency noise. The filter cutoff frequency can be set to 100 Hz to 200 Hz. Then, rising edge detection is performed on the filtered pressure signal: the pressure value at the current sampling point is compared with the pressure value at the previous sampling point. If the current value is greater than the previous value and the difference exceeds a preset noise threshold, it is recorded as a rising edge. The noise threshold can be set according to the background noise level of the pressure signal, generally between 0.5% and 1.0% of the pressure sensor's range. The total number of rising edges is counted per unit time, and this number is recorded as the pressure impact frequency, measured in times per minute. The pressure impact frequency reflects the frequency of pressure impacts on the propulsion hydraulic system: during hard rock excavation, the cutter head generates frequent impact loads, increasing the number of rising edges in the hydraulic system pressure signal; during soft rock excavation, the load is stable, and the number of rising edges in the pressure signal decreases.

[0088] The second harmonic ratio (HHR) is multiplied by the pressure impact frequency, and the product is recorded as the geological feature fusion value. The HHR is dimensionless, while the pressure impact frequency is measured in cycles per minute (CPM); therefore, the geological feature fusion value is also measured in CPM. This fusion value integrates information from both the electrical and hydraulic domains: the HHR reflects the degree of nonlinear distortion of the rock-breaking load waveform in the frequency domain, while the pressure impact frequency reflects the impact frequency of the rock-breaking load in the time domain. Multiplying the two helps enhance the signal characteristics of lithological changes, thus improving the reliability of lithological identification.

[0089] One rotation of the cutterhead constitutes one acquisition cycle. At the end of each acquisition cycle, the geological feature fusion value for that cycle is calculated. The geological feature fusion value of the current acquisition cycle is compared with the geological feature fusion values ​​of the previous two acquisition cycles to determine the trend of this value over three consecutive acquisition cycles. If the geological feature fusion value shows a monotonically increasing trend over three consecutive acquisition cycles (i.e., the value in the third cycle is greater than the value in the second cycle, and the value in the second cycle is greater than the value in the first cycle), it is determined that the lithology of the tunnel face is hardening. In this case, the average thrust extracted from the historical tunneling dataset in step two is multiplied by a first correction coefficient greater than 1, and the product is used as the corrected average thrust input to the Long Short-Term Memory network model. The first correction coefficient can be set to a value between 1.05 and 1.30; for example, it can be set to 1.10 when homogeneous strata transition to hard rock, and to 1.20 when there is a sudden change in lithology. If the geological feature fusion value shows a monotonically decreasing trend over three consecutive acquisition cycles (i.e., the value in the third cycle is less than the value in the second cycle, and the value in the second cycle is less than the value in the first cycle), it is determined that the lithology of the tunnel face is softening. At this point, the mean thrust extracted in step two is multiplied by a second correction factor less than 1, and the product is used as the corrected mean thrust input to the long short-term memory network model. The second correction factor can be set to a value between 0.70 and 0.95; for example, it can be set to 0.90 when transitioning from homogeneous strata to soft rock, and to 0.80 when there is a sudden change in lithology. If the geological feature fusion value does not show a monotonically increasing or decreasing trend over three consecutive acquisition cycles, it is determined that the lithology has not changed significantly, and no correction is made to the mean thrust; the original mean thrust is directly input into the long short-term memory network model.

[0090] Compared to existing technologies that rely primarily on cutterhead vibration signals or manual observation of rock debris for geological assessment, this invention obtains the second harmonic ratio by collecting the current signal from the cutterhead drive motor and the pressure impact frequency by collecting the pressure signal from the propulsion hydraulic system. Multiplying the two yields a fused geological feature value, which is then used to correct the average input propulsion force of the long short-term memory network model. This allows the acquisition of geological information to be independent of additional vibration sensors or manual observation. Furthermore, both the current and pressure signals are inherent and usable signals during TBM tunneling, eliminating the need for additional dedicated geological exploration equipment.

[0091] According to another embodiment of the present invention, the intelligent control method for TBM construction specifies the acquisition of current and pressure signals, the construction of sliding windows, and the calculation of geological feature fusion values. The current signal of the cutterhead drive motor is acquired through a current transformer or a Hall effect current sensor installed on the stator winding lead of the motor. The current sensor can be a closed-loop Hall effect current sensor, with its measurement range selected according to the rated current of the motor, generally 1.2 to 2.0 times the rated current, and the output signal is a 4 mA to 20 mA current signal or a 0 V to 10 V voltage signal. The pressure signal of the propulsion hydraulic system is acquired through a pressure sensor installed near the oil inlet of the propulsion hydraulic cylinder. The pressure sensor can be a piezoresistive or strain gauge pressure sensor, with its range selected according to the rated working pressure of the propulsion hydraulic system, generally 1.5 to 2.0 times the rated pressure, and the output signal is also a 4 mA to 20 mA current signal or a 0 V to 10 V voltage signal. The output signals of the two sensors are filtered and amplified by their respective signal conditioning circuits before being sent to the multi-channel data acquisition card built into the industrial control computer. The sampling frequency of the data acquisition card can be set from 1 kHz to 10 kHz to ensure that the components of the cutter head rotation frequency and its integer multiples in the current amplitude envelope signal, as well as the transient characteristics of the rising edge in the pressure signal, can be completely preserved.

[0092] One acquisition cycle is defined as one revolution of the cutter head. The duration of one revolution is determined by the current cutter head rotation speed. For example, when the cutter head rotation speed is 4 revolutions per minute, one acquisition cycle lasts approximately 15 seconds; when the cutter head rotation speed is 6 revolutions per minute, one acquisition cycle lasts approximately 10 seconds. Within each acquisition cycle, the data acquisition card continuously acquires current and pressure signals, and the acquired data is temporarily stored in the buffer of the industrial control computer's memory in chronological order.

[0093] At the end of each acquisition cycle, the current signal from the current acquisition cycle and the four preceding consecutive acquisition cycles (a total of five acquisition cycles) is retrieved from the buffer and stored in a current sliding window. The current sliding window is a data storage area allocated in the industrial control computer's memory. Its capacity varies with the total number of sampling points in the five acquisition cycles at the current cutter head speed, meaning it always accommodates all current sampling data from the current acquisition cycle and the four preceding consecutive acquisition cycles. Five acquisition cycles are chosen because the frequency resolution of the Fourier transform depends on the duration of the signal involved in the transform; the signal length of five acquisition cycles is sufficient to obtain an envelope spectrum with adequate frequency resolution, effectively distinguishing the fundamental frequency component and its second harmonic component of the cutter head rotation frequency. Simultaneously, the pressure signal from the current acquisition cycle and the two preceding consecutive acquisition cycles (a total of three acquisition cycles) is retrieved from the buffer and stored in a pressure sliding window. The pressure sliding window is also a data storage area allocated in the industrial control computer's memory, and its capacity varies with the total number of sampling points in the three acquisition cycles at the current cutter head speed, meaning it always accommodates all pressure sampling data from the current acquisition cycle and the two preceding consecutive acquisition cycles. The frequency of pressure shocks is a statistical feature in the time domain. If the window is too long, it will smooth out the transient shock information, and if the window is too short, the statistical sample will be insufficient. Three acquisition cycles can achieve a reasonable balance between statistical sufficiency and transient preservation.

[0094] Both the current sliding window and the pressure sliding window use a first-in, first-out (FIFO) update method. At the end of each new acquisition cycle, the data for the current acquisition cycle is stored at the end of the window, while the data from the oldest acquisition cycle at the beginning of the window is removed, thus ensuring that the window always contains data from the latest five or three acquisition cycles.

[0095] At the end of each acquisition cycle, the current signal from all five acquisition cycles within the current sliding window is envelope-detected and then subjected to Fast Fourier Transform (FFT) to obtain the envelope spectrum. Before performing the FFT, the current amplitude envelope signal can be de-biased and a Hanning or Hamming window can be applied to reduce spectral leakage. The number of points in the FFT should not be less than the signal length within the current sliding window; generally, a number of points greater than or equal to the signal length and a power of 2 is selected. For example, when the sampling frequency is 2 kHz and the cutter head speed is 4 revolutions per minute, the signal length within the current sliding window is approximately 150,000 sampling points, in which case 262,144 points can be selected. In the envelope spectrum obtained after the transformation, the horizontal axis represents frequency, and the vertical axis represents amplitude, with the amplitude unit being amperes.

[0096] Extract the fundamental frequency and second harmonic amplitude of the cutter head rotation frequency from the envelope spectrum. The fundamental frequency of the cutter head rotation frequency is equal to the cutter head rotation speed divided by 60, and the unit is Hertz (Hz). For example, when the cutter head rotation speed is 4 revolutions per minute, the fundamental frequency is approximately 0.067 Hz, and the second harmonic is approximately 0.133 Hz. Since the envelope spectrum is discrete, the fundamental frequency and second harmonic may not fall exactly at discrete frequency points. In this case, the local maximum amplitude can be searched within the frequency range near the fundamental frequency and the second harmonic, and this maximum value can be taken as the amplitude of that frequency component. Calculate the ratio of the second harmonic amplitude to the fundamental frequency amplitude. This ratio is dimensionless and is denoted as the second harmonic ratio of the current sliding window.

[0097] Rising edge statistics are performed on the pressure signal across all three acquisition cycles within the pressure sliding window. First, the pressure signal is low-pass filtered to eliminate high-frequency noise; the filter cutoff frequency can be set to 100 Hz to 200 Hz. Then, rising edge detection is performed: the pressure value at the current sampling point is compared with the pressure value at the previous sampling point. If the current value is greater than the previous value and the difference exceeds a preset noise threshold, it is recorded as a rising edge. The noise threshold can be set according to the background noise level of the pressure signal, typically between 0.5% and 1.0% of the pressure sensor's range. The total number of rising edges within the pressure sliding window is counted, and this number is divided by the total duration of the three acquisition cycles (in minutes) to obtain the number of rising edges per unit time, recorded as the pressure impact frequency of the pressure sliding window, in times per minute.

[0098] The second harmonic ratio of the current sliding window is multiplied by the pressure impact frequency of the pressure sliding window, and the product is used as the geological feature fusion value for that acquisition period. The second harmonic ratio is dimensionless, and the pressure impact frequency is measured in cycles per minute; therefore, the geological feature fusion value is also measured in cycles per minute. The current sliding window covers five acquisition periods, ensuring the frequency resolution required for envelope spectrum analysis; the pressure sliding window covers three acquisition periods, ensuring the sensitivity of pressure impact frequency statistics to transient information. The end times of both windows are aligned with the end of the current acquisition period; therefore, the two multiplied feature values ​​reflect the rock-breaking process within the same time period.

[0099] Compared to existing technologies that use the same time window for feature extraction of current and pressure signals, this invention sets sliding windows of different lengths to suit the frequency domain analysis characteristics of current signals and the time domain statistical characteristics of pressure signals. The current window has five acquisition cycles and the pressure window has three acquisition cycles, so that the length of the two windows matches the time scale requirements of their respective feature extraction methods, thereby achieving better results in both frequency resolution and transient information preservation.

[0100] According to another embodiment of the present invention, the intelligent control method for TBM construction further defines the alignment of the current sliding window and the pressure sliding window on the time axis, as well as the specific calculation process of the geological feature fusion value. At the end of each acquisition cycle, the current sliding window contains the current signal of the current acquisition cycle and the previous four consecutive acquisition cycles, totaling five acquisition cycles, with its time axis coverage ranging from the start time of the earliest acquisition cycle to the end time of the current acquisition cycle. The pressure sliding window contains the pressure signal of the current acquisition cycle and the previous two consecutive acquisition cycles, totaling three acquisition cycles, with its time axis coverage ranging from the start time of the earliest acquisition cycle to the end time of the current acquisition cycle. Although the two windows have the same end time, both aligned with the end time of the current acquisition cycle, their start times are different; the start time of the current sliding window is two acquisition cycles earlier than the start time of the pressure sliding window. Therefore, the intersection of the two sliding windows on the time axis is the three acquisition cycles covered by the pressure sliding window, namely the current acquisition cycle and the previous two consecutive acquisition cycles. This intersection is defined as a common time window, the length of which is equal to the length of the pressure sliding window, i.e., three acquisition cycles.

[0101] At the end of each acquisition cycle, envelope detection is performed on the current signal from all five acquisition cycles within the current sliding window, followed by Fast Fourier Transform (FFT) to obtain the envelope spectrum. Before the FFT, the current amplitude envelope signal can be de-DC biased and a Hanning or Hamming window can be applied to reduce spectral leakage. The fundamental frequency and second harmonic amplitude of the cutter head rotation frequency are extracted from the envelope spectrum, and the ratio of the second harmonic amplitude to the fundamental frequency amplitude is calculated. This ratio is denoted as the second harmonic ratio of the current sliding window. The calculation of the second harmonic ratio is still based on the current signal from all five acquisition cycles within the current sliding window, because the frequency resolution of the Fourier transform depends on the duration of the signal involved in the transform, and the signal length of five acquisition cycles provides sufficient frequency resolution to effectively distinguish the fundamental and second harmonic components. Since the acquisition period is inversely proportional to the cutter head rotation speed, the lower the cutter head rotation speed, the longer the single acquisition period, the total duration of the current sliding window increases accordingly, and the frequency resolution is improved accordingly. The relative interval between the fundamental frequency and the second harmonic in the envelope spectrum remains stable. Therefore, this feature has a consistent frequency resolution capability within the range of cutter head rotation speed variation.

[0102] For pressure signals within a common time window (i.e., the pressure signals of the current acquisition cycle and the two preceding consecutive acquisition cycles, totaling three acquisition cycles), rising edge statistics are performed. First, the pressure signals within the common time window are low-pass filtered to eliminate high-frequency noise; the filter cutoff frequency can be set to 100 Hz to 200 Hz. Then, rising edge detection is performed: the pressure value at the current sampling point is compared with the pressure value at the previous sampling point. If the current value is greater than the previous value and the difference exceeds a preset noise threshold, it is recorded as a rising edge. The noise threshold can be set according to the background noise level of the pressure signal, typically between 0.5% and 1.0% of the pressure sensor's range. The total number of rising edges within the common time window is counted, and this number is divided by the total duration of the common time window (in minutes) to obtain the number of rising edges per unit time, recorded as the pressure impact frequency of the common time window, in times per minute. Since the common time window only covers three acquisition cycles, the pressure impact frequency statistics are strictly based on the pressure signals within this time period and do not involve pressure data from the earliest two acquisition cycles outside the current sliding window.

[0103] The second harmonic ratio of the current sliding window is multiplied by the pressure impact frequency of the common time window, and the product is used as the geological feature fusion value for that acquisition cycle. The second harmonic ratio is calculated based on the current signal of five acquisition cycles, ensuring the frequency resolution required for spectral analysis; the pressure impact frequency is calculated based on the pressure signal of three acquisition cycles, avoiding the smoothing out of transient impact information due to an excessively long window. When multiplied, the second harmonic ratio reflects the average load nonlinearity over the five acquisition cycles, while the pressure impact frequency reflects the average impact frequency over the three acquisition cycles. Since the three acquisition cycles of the common time window are completely contained within the five acquisition cycles of the current sliding window, and the end times of both windows are aligned with the end of the current acquisition cycle, the two multiplied feature values ​​have a clear overlap on the time axis, avoiding geological information aliasing caused by incomplete overlap of time periods.

[0104] Compared to existing technologies that directly fuse signals from different sources after feature extraction using the same time window, this invention determines the intersection of the current sliding window and the pressure sliding window as a common time window. This ensures that the statistics of pressure impact frequency are strictly based on the time period covered by both windows, eliminating the interference of the two additional acquisition cycles in the current sliding window on time alignment. This results in a precise correspondence between the two multiplied feature values ​​on the time axis, which is beneficial to improving the accuracy of geological feature fusion values.

[0105] According to another embodiment of the present invention, the specific method for outputting the thrust adjustment amount and torque adjustment amount of the long short-term memory network model in step two is refined in the provided intelligent control method for TBM construction.

[0106] The Long Short-Term Memory (LSTM) network model is deployed on the processor of the TBM industrial computer. At each time step, after receiving input features, the model performs forward computation and outputs a predicted thrust and a predicted torque value from its output layer. The unit of the predicted thrust is kilonewtons (kN), and the unit of the predicted torque is kilonewton-meters (kNm). Each time step corresponds to a sampling period, which can be set to 1 second or any value between 0.5 and 2 seconds depending on actual needs. At each time step, the industrial computer simultaneously acquires the actual thrust and torque values. The actual thrust is obtained by converting the pressure sensor installed at the inlet of the hydraulic cylinder into the effective working area of ​​the hydraulic cylinder. The actual torque is calculated from the current and speed of the cutter head drive motor or directly measured by a torque sensor installed on the cutter head drive shaft.

[0107] The instantaneous thrust error is obtained by subtracting the predicted thrust value from the actual thrust value collected at the current time step. Similarly, the instantaneous torque error is obtained by subtracting the predicted torque value from the actual torque value collected at the current time step. The unit for the instantaneous thrust error is kilonewtons (kN), and the unit for the instantaneous torque error is kilonewton-meters (kNm). A positive instantaneous error indicates that the model's predicted value is greater than the actual value, while a negative instantaneous error indicates that the model's predicted value is less than the actual value.

[0108] The instantaneous thrust errors of P consecutive time steps, including the current time step, are arranged in chronological order to form a thrust error sequence. P is a preset positive integer, which can be selected between 10 and 50, for example, 20. The value of P determines the length of the error sequence used to determine the prediction delay; a P value that is too small results in insufficient statistical samples, while a P value that is too large reduces the sensitivity to recent changes. The variance of all elements in the thrust error sequence is calculated and denoted as the thrust error variance. The thrust error variance is calculated as follows: first, the arithmetic mean of the thrust error sequence is calculated; then, the square of the difference between each element and this mean is calculated; the sum of all squared values ​​is divided by the sequence length, which is the thrust error variance σ. 2 The calculation method is as follows: , where e i Let be the error value at the i-th time step in the propulsion error sequence. Let P be the arithmetic mean of the propulsion error sequence, and P be the length of the error sequence. The variance reflects the degree of dispersion of the error sequence around its mean. Similarly, the torque error sequence is constructed by taking the instantaneous torque errors of P consecutive time steps, including the current time step, and the variance of the torque error sequence is calculated and denoted as the torque error variance.

[0109] When the thrust error variance exceeds the first variance threshold, the current thrust prediction error is determined to be outside the allowable range, requiring correction. In this case, the thrust prediction value at the current time step is replaced with the arithmetic mean or weighted moving average of the actual thrust values ​​from a preset number of time steps prior to the current time step. The replaced value is output as the thrust adjustment amount. This preset number is a preset positive integer less than P, and its specific value is determined based on the estimated delay steps: 1 time step is used when the estimated delay steps are 0 or 1; the estimated delay steps are used when the estimated delay steps are greater than 1 and less than or equal to the first delay threshold; and a preset fixed number of time steps are used when the estimated delay steps are greater than the first delay threshold. This preset number being less than P ensures that the historical actual value window used to replace the prediction value is shorter than the error sequence window used to determine the delay, ensuring that the replacement value is based on updated actual data. The arithmetic mean of the actual thrust values ​​is calculated by summing the actual thrust values ​​from each of the preset number of time steps prior to the current time step and dividing by this preset number. When the thrust error variance is less than or equal to the first variance threshold, the current thrust prediction error is determined to be within the allowable range, and the thrust prediction value of the current time step is directly output as the thrust adjustment amount.

[0110] Similarly, when the torque error variance exceeds the second variance threshold, the current torque prediction error is determined to be outside the allowable range, and the prediction value needs to be corrected. The torque prediction value at the current time step is replaced with the arithmetic mean or weighted moving average of the actual torque values ​​at the previous preset number of time steps, and the replaced value is output as the torque adjustment amount. When the torque error variance is less than or equal to the second variance threshold, the torque prediction value at the current time step is directly output as the torque adjustment amount.

[0111] The first and second variance thresholds can be preset based on the statistical distribution of propulsion and torque prediction errors under normal operating conditions. The first variance threshold can be set within the range of 2 to 5 times the historical propulsion error variance, for example, it can be set to 3 times the historical propulsion error variance. The second variance threshold can be set within the range of 2 to 5 times the historical torque error variance, for example, it can be set to 3 times the historical torque error variance. The historical propulsion error variance and torque error variance can be calculated by collecting prediction error data for a period of time during the initial stable tunneling phase after TBM startup.

[0112] Compared to existing technologies that directly use the original predicted values ​​of the Long Short-Term Memory network model as control commands, this invention determines the reliability of the current predicted value by calculating the variance of the prediction error sequence. When the error variance exceeds a threshold, the current predicted value is replaced by the arithmetic mean of the actual values ​​over the most recent period. This allows the control commands to still be generated based on the latest actual tunneling status when the model generates delayed predictions, thereby reducing the impact of prediction delay on the control effect.

[0113] According to another embodiment of the present invention, the intelligent control method for TBM construction further refines the specific method of replacing the predicted value with the historical actual value when the variance of the thrust error is greater than the first variance threshold. When the variance of the thrust error is greater than the first variance threshold, the variance of the thrust error sequence for P consecutive time steps, including the current time step, is first calculated and denoted as the current variance. The calculation method of the current variance is the same as the aforementioned calculation method of the thrust error variance, that is, first the arithmetic mean of the thrust error sequence is calculated, then the square of the difference between each error value and the mean is calculated, and all squared values ​​are added together and divided by the sequence length.

[0114] Subsequently, starting from k=1, the values ​​are sequentially selected up to K, where K is the preset maximum delay step. The value of K can be determined based on the maximum possible response delay time of the propulsion hydraulic system, and is generally set to an integer between 5 and 20, for example, 10. For each k value, the difference between the predicted propulsion force value of P consecutive time steps including the current time step and the actual propulsion force value offset k time steps before the current time step is calculated. Specifically, when k=1, the predicted propulsion force value of P consecutive time steps including the current time step is subtracted from the actual propulsion force value of the previous time step, resulting in P differences; when k=2, the predicted propulsion force value of these P time steps is subtracted from the actual propulsion force value of the previous two time steps, resulting in P differences; and so on, until k=K. The P differences corresponding to each k value are arranged in chronological order to form the offset error sequence corresponding to that k value. The length of this offset error sequence is the same as the length of the propulsion error sequence, which is P time steps. Calculate the variance of each offset error sequence separately, denoted as the offset variance. The method for calculating the offset variance is the same as that for calculating the current variance. Furthermore, k takes the value of an integer from 0 to K; when k equals 0, the offset error sequence is the thrust error sequence, and its offset variance is equal to the current variance.

[0115] The current variance is compared with each offset variance one by one. When there is an offset variance smaller than the current variance, the k value corresponding to the smallest offset variance is selected as the estimated delay step number. When all offset variances are greater than or equal to the current variance, the estimated delay step number is 0. Determined in the following ways: ,in, To estimate the number of delay steps, The current variance, Let k be the offset variance with an offset of k time steps. The offset variance with an offset of 0 time steps is equal to the current variance. The value of k that minimizes the offset variance is selected as the estimated delay step number. The principle is as follows: when the model prediction produces a delay, the fluctuation of the error between the predicted value and the actual value after offsetting the predicted value by k time steps is less than the fluctuation when there is no offset. At this time, the offset variance corresponding to this offset is the minimum. When the prediction does not produce a delay, any offset will increase the fluctuation of the error, and all offset variances are not less than the current variance, so the estimated delay step number is 0. By traversing all possible values ​​of k and finding the k with the smallest offset variance, the current prediction delay step number can be estimated.

[0116] After obtaining the estimated delay steps, different processing methods are applied based on the specific value of the estimated delay steps. When the estimated delay steps are 0 or 1, it indicates that the current prediction delay does not exceed one sampling period. In this case, it is counted as one sampling period, and the predicted thrust value of the current time step is replaced with the actual thrust value of the previous time step, and the result is used as the thrust adjustment amount output. When the estimated delay steps are greater than 1 and less than or equal to a preset first delay threshold, it indicates that a slight delay has occurred in the prediction. The predicted thrust value of the current time step is replaced with the arithmetic mean of the actual thrust values ​​of the estimated delay steps preceding the current time step, and the replaced value is used as the thrust adjustment amount output. The first delay threshold can be set to an integer between 3 and 6, for example, it can be set to 3, and the first delay threshold is less than K. When the estimated delay steps are greater than the first delay threshold, it indicates that a severe delay has occurred in the prediction. The predicted thrust value of the current time step is replaced with the weighted moving average of the actual thrust values ​​of the preset fixed number of time steps preceding the current time step, and the replaced value is used as the thrust adjustment amount output. The preset fixed number of time steps can be selected between 5 and 15, for example, it can be set to 10, and the preset fixed number is less than P. The weighted moving average refers to the sum of the actual values ​​of each propulsion force within the preset fixed number of time steps preceding the current time step, assigning different weights to each. Time steps closer to the current time step are assigned higher weights, and time steps farther away are assigned lower weights, with the sum of all weights equal to 1. That is, the weighted moving average. The calculation method is as follows: Where M is a preset fixed number of time steps, x i w represents the actual thrust value at the i-th time step prior to the current time step (i=1 being the most recent time step). i The weights are assigned accordingly, with time steps closer to the current time step receiving a larger weight and time steps farther away receiving a smaller weight. For example, a linearly decreasing weighting can be used, with the time step closest to the current time step having the largest weight and the farthest time step having the smallest weight.

[0117] Compared to existing technologies that only determine whether a delay has occurred and uniformly apply the same compensation method, this invention estimates the specific number of delay steps by comparing the current variance with each offset variance, and adopts graded processing according to the different delay steps. For mild delays, the arithmetic mean is used instead, and for severe delays, the weighted moving average is used instead, so that the compensation method matches the degree of delay.

[0118] According to another embodiment of the present invention, the provided intelligent control method for TBM construction further refines the specific methods for calculating the variance of the propulsion error sequence over P consecutive time steps, including the current time step, and for calculating the variance of each offset error sequence. When a trend component exists in the propulsion error sequence, for example, when the propulsion demand gradually increases due to tool wear and the model fails to track this trend in time, the mean of the error sequence will drift slowly over time. In this case, even if the prediction is not delayed, the variance of the error sequence will increase due to the influence of the trend component, thereby interfering with the estimation result of the delay step number. To eliminate this interference, the error sequence is first subjected to linear detrending processing before calculating the variance.

[0119] For the propulsion error sequence E=[e1,e2,…,e P ](e i To perform linear regression on the error value at the i-th time step, we treat each element in the thrust error sequence as the dependent variable and the corresponding time step number i as the independent variable, and fit a straight line. The slope *b* of the straight line is obtained. This slope reflects the overall trend direction and rate of change of the error sequence over time. A positive slope indicates an upward trend in the error sequence, while a negative slope indicates a downward trend. A trend term sequence T = [a+b×1, a+b×2, ..., a+b×P] is generated based on the slope *b* and intercept *a* of the linear regression. Each element in the trend term sequence is equal to the sum of the intercept *a* and the slope *b* multiplied by the corresponding time step number. The detrended propulsion error sequence is obtained by subtracting the trend value at the corresponding time from each value in the propulsion error sequence. Detrending eliminates the linear trend component in the error sequence, causing the detrended error sequence to fluctuate around zero.

[0120] Calculate the variance of the detrended driving force error sequence as the current variance. Similarly, perform the same linear detrending process on each offset error sequence, i.e., perform linear regression on each offset error sequence to obtain its own slope and intercept, generate its own trend term sequence based on its slope and intercept, subtract the trend value of the corresponding time in the respective trend term sequence from each value of each offset error sequence to obtain each detrended offset error sequence, and then calculate the variance of each detrended offset error sequence as the offset variance. Compare the current variance after detrending with each offset variance after detrending one by one. When there is an offset variance smaller than the current variance, select the k value corresponding to the smallest offset variance as the estimated delay step; when all offset variances are greater than or equal to the current variance, the estimated delay step is 0.

[0121] Detrending processing ensures that variance calculation truly reflects the degree of fluctuation of the error sequence around its mean, rather than the overall trend change of the error sequence, thus improving the accuracy of delay step estimation. The setting methods for parameters such as the first variance threshold, the second variance threshold, and the first delay threshold are the same as those described above and will not be repeated here.

[0122] Compared to existing technologies that directly calculate the variance of the error sequence, this invention performs linear detrending processing on the error sequence before calculating the variance. This eliminates the interference of trend components caused by factors such as tool wear and slow changes in the surrounding rock on the variance calculation, allowing the variance to more accurately reflect the degree of fluctuation of the prediction error around the mean, thereby improving the accuracy of the delay step estimation.

[0123] According to another embodiment of the present invention, in the provided intelligent control method for TBM construction, the Long Short-Term Memory (LSTM) network model is deployed on the processor of the TBM industrial control computer and is updated online during the tunneling process in the following manner: An update cycle is defined as N rotations of the cutterhead, where N is an integer greater than 1. The value of N can be selected between 5 and 20, for example, it can be set to 10. The number of cutterhead rotations is counted by a rotary encoder or proximity switch installed on the cutterhead drive shaft. When the cumulative number of rotations reaches N, a model update judgment is triggered. The update cycle is set in accordance with the physical tunneling progress of the TBM, rather than with a fixed time period. When the TBM tunnels at a faster speed in homogeneous strata, the time required for the cutterhead to rotate N times is shorter, and the model update frequency is correspondingly increased; when the tunneling speed is slower in hard rock strata, the time required for the cutterhead to rotate N times is longer, and the model update frequency is correspondingly decreased. This binding method allows the model update rhythm to naturally synchronize with the changing rhythm of the tunneling conditions.

[0124] At the end of each update cycle, the predicted thrust and torque values ​​output by the Long Short-Term Memory (LSTM) network model at each time step within that update cycle are read from the industrial control computer's memory. Simultaneously, the actual thrust and torque value sequences corresponding to each time step are read from the data acquisition module. The actual thrust value is obtained by converting the effective working area of ​​the hydraulic cylinder using a pressure sensor installed at the inlet of the thrust hydraulic cylinder. The actual torque value is calculated using the current and speed of the cutter head drive motor or directly measured by a torque sensor installed on the cutter head drive shaft.

[0125] After obtaining the thrust prediction value sequence, the difference between the predicted values ​​at two adjacent time steps is calculated. Specifically, the predicted value at the second time step is subtracted from the predicted value at the first time step, the predicted value at the third time step is subtracted from the predicted value at the second time step, and so on, until the end of the sequence, resulting in the thrust prediction difference sequence. Each element in the difference sequence reflects the direction of change of the predicted value between two adjacent time steps: a positive difference indicates an increase in the predicted value, a negative difference indicates a decrease in the predicted value, and a difference of zero indicates that the predicted value remains unchanged. Similarly, the difference between the actual values ​​at two adjacent time steps in the actual thrust value sequence is calculated to obtain the actual thrust difference sequence.

[0126] The thrust prediction difference sequence is compared element by element at corresponding positions with the actual thrust difference sequence, and the number of elements with the same sign is counted. Same sign includes three cases: both are positive, both are negative, or both are zero. The number of elements with the same sign is divided by the total length of the difference sequence to obtain the thrust trend consistency rate. The trend consistency rate ranges from 0 to 1, or can be expressed as 0% to 100%. The trend consistency rate reflects the extent to which the direction of change of the model's predicted value is consistent with the direction of change of the actual value. When the model prediction is accurate, the direction of change of the predicted value should be consistent with the direction of change of the actual value, resulting in a high trend consistency rate. When the model introduces a prediction delay, the direction of change of the predicted value lags behind the direction of change of the actual value; for example, the actual value has started to rise while the predicted value is still falling. In this case, the difference signs are opposite, and the trend consistency rate decreases.

[0127] Similarly, the difference between two adjacent time-step predicted values ​​in the torque prediction value sequence is calculated to obtain the torque prediction difference sequence. The difference between two adjacent time-step actual values ​​in the actual torque value sequence is calculated to obtain the actual torque difference sequence. The sign consistency rate between the torque prediction difference sequence and the actual torque difference sequence is calculated and denoted as the torque trend consistency rate. The calculation method is the same as the thrust trend consistency rate, that is, the number of elements with the same sign at corresponding positions in the torque prediction difference sequence and the actual torque difference sequence is counted and divided by the total length of the difference sequences.

[0128] The thrust trend consistency rate is compared with a preset first consistency rate threshold, and the torque trend consistency rate is compared with a preset second consistency rate threshold. The first and second consistency rate thresholds can be preset based on the model's prediction performance under normal operating conditions; for example, both can be set to values ​​between 0.70 and 0.90, specifically 0.80. When the thrust trend consistency rate is lower than the first consistency rate threshold, it indicates that the direction of change in the predicted thrust value has deviated significantly from the direction of change in the actual value, and the model's prediction of thrust can no longer accurately reflect the actual trend, requiring a model update. Similarly, when the torque trend consistency rate is lower than the second consistency rate threshold, it also indicates that the torque prediction can no longer accurately reflect the actual trend, requiring a model update. When either the thrust trend consistency rate is lower than the first or the torque trend consistency rate is lower than the second consistency rate threshold, all data frames within that update period are added as new training samples to the training set of the Long Short-Term Memory (LSTM) network model. The LTM network model is then incrementally trained, and the model parameters are updated. The incremental training method is as follows: New training samples are merged with the existing training set, and several rounds of backpropagation training are performed based on the original model parameters to adapt the model parameters to the new data distribution. A lower learning rate is used to avoid catastrophic forgetting. When the consistency rate of the propulsion trend is greater than or equal to the first consistency rate threshold and the consistency rate of the torque trend is greater than or equal to the second consistency rate threshold, it indicates that the model's current prediction trend is still accurate, and no update is triggered, keeping the model parameters unchanged.

[0129] Compared to existing technologies that trigger model updates based on the absolute value of the prediction error or a fixed period, this invention determines the accuracy of the model's predicted trend by calculating the sign consistency rate between the predicted difference sequence and the actual difference sequence. The trend consistency rate only focuses on whether the direction of change of the predicted value is consistent with the actual value. It is not sensitive to the dimensions and magnitude, and is not affected by changes in the magnitude of propulsion force and torque. It can trigger updates at an early stage when the model's predicted trend begins to deviate from the actual trend, thereby improving the model's timeliness in responding to changes in operating conditions.

[0130] According to another embodiment of the present invention, the intelligent control method for TBM construction further refines the specific method for determining whether the model has degraded. When the thrust trend consistency rate is greater than or equal to the first consistency rate threshold, it indicates that the overall direction of change of the thrust prediction value is consistent with the actual value, and the model's direction prediction is accurate. However, accurate direction does not guarantee accurate amplitude. The model's predicted value may have a continuous amplitude deviation from the actual value. For example, the predicted value is always 10% higher than the actual value, and this deviation gradually increases over time. This amplitude drift is also a manifestation of model degradation and requires an update. To detect this amplitude drift, the relative deviation between the thrust prediction value and the actual thrust value at each time step within the update cycle is further calculated. The relative deviation is calculated by subtracting the actual thrust value from the thrust prediction value and then dividing by a non-zero actual thrust value. The relative deviation is dimensionless. The relative deviations at each time step are arranged in chronological order to obtain a thrust relative deviation sequence.

[0131] Calculate the slope of the linear regression for the relative deviation sequence of propulsion force. The specific method of linear regression is to treat each element in the relative deviation sequence as the dependent variable and the corresponding time sequence number as the independent variable, fitting a straight line. The slope of this line is the linear regression slope. The unit of the linear regression slope is per time step, reflecting the overall trend of the relative deviation throughout the entire update cycle. A positive slope indicates that the relative deviation is increasing over time, meaning the deviation between the model's predicted value and the actual value is continuously increasing; a negative slope indicates that the relative deviation is decreasing over time, meaning the deviation between the model's predicted value and the actual value is continuously decreasing.

[0132] The linear regression slope is compared with preset positive slope thresholds and preset negative slope thresholds. The preset positive slope threshold is a positive value, which can be set between 0.01 and 0.05, for example, 0.02. The preset negative slope threshold is a negative value, which can be set between -0.05 and -0.01, for example, -0.02. When the linear regression slope is positive and greater than the preset positive slope threshold, it is determined that the predicted propulsion value has a positive amplitude drift relative to the actual value, meaning the predicted value is continuously larger than the actual value. In this case, all data frames within the update period are added as new training samples to the training set of the Long Short-Term Memory (LSTM) network model, and the model parameters are updated after incremental training. When the linear regression slope is negative and its absolute value is greater than the absolute value of the preset negative slope threshold, it is determined that a negative amplitude drift has occurred, meaning the predicted value is continuously smaller than the actual value, and incremental training is also triggered. When the absolute value of the linear regression slope is less than or equal to the corresponding threshold, it is determined that no significant amplitude drift has occurred, and no update is triggered.

[0133] Compared to existing technologies that rely solely on trend consistency rate to determine whether a model needs updating, this invention further calculates the linear regression slope of the relative deviation sequence, provided that the trend consistency rate meets the standard. By determining the sign and magnitude of the slope, it judges whether there is a persistent amplitude drift between the predicted value and the actual value, enabling the model update to respond to both directional deviation and amplitude drift degradation modes simultaneously.

[0134] According to another embodiment of the present invention, the intelligent control method for TBM construction further refines the specific method for detecting abrupt amplitude drift. This method complements the aforementioned linear regression slope method. The relative deviation sequence of propulsion force is divided into a first segment sequence and a second segment sequence in chronological order. The division method is as follows: if the length of the relative deviation sequence is even, the first segment sequence and the second segment sequence each occupy half; if the length is odd, the middle element can be assigned to either the first segment sequence or the second segment sequence, for example, it can be assigned to the first segment sequence. The arithmetic mean of the first segment sequence is calculated as the first segment mean, and the arithmetic mean of the second segment sequence is calculated as the second segment mean. Both the first segment mean and the second segment mean are dimensionless, reflecting the average level of the relative deviation in the first and second halves of the update cycle.

[0135] Calculate the difference between the mean of the latter half and the mean of the former half. A positive difference indicates that the average level of the relative deviation in the latter half is higher than that in the first half, meaning the deviation of the predicted value from the actual value increases in the latter half of the update cycle; a negative difference indicates that the average level of the relative deviation in the latter half is lower than that in the first half, meaning the deviation decreases in the latter half. This difference reflects the magnitude and direction of the change in relative deviation from the first half to the second half within the update cycle and is highly sensitive to abrupt amplitude drift. When amplitude drift occurs near the middle of the update cycle, the difference between the mean of the former half and the mean of the latter half can be immediately reflected.

[0136] The difference between the mean of the later segment and the mean of the earlier segment is compared with two preset difference thresholds: a first difference threshold and a second difference threshold. The first difference threshold is a positive value, which can be set between 0.05 and 0.15, for example, 0.10. The second difference threshold is a negative value, which can be set between -0.15 and -0.05, for example, -0.10. When the difference is greater than the first difference threshold, a positive abrupt amplitude shift is determined to have occurred, meaning that the deviation of the predicted value from the actual value increases significantly in the latter half of the update cycle. All data frames within that update cycle are added as new training samples to the training set of the Long Short-Term Memory (LSTM) network model, and the model parameters are updated after incremental training. When the difference is less than the second difference threshold, a negative abrupt amplitude shift is determined to have occurred, and incremental training is also triggered. When the difference is greater than or equal to the second difference threshold and less than or equal to the first difference threshold, no significant abrupt amplitude shift is determined to have occurred, and no update is triggered.

[0137] The piecewise mean comparison method detects the difference in mean between two segments before and after the update cycle, and is sensitive to abrupt amplitude drift; while the linear regression slope method detects the overall trend of change throughout the entire update cycle, and is sensitive to gradual amplitude drift. The two methods detect amplitude drift from different perspectives, thus complementing each other. In practical applications, both methods can be used simultaneously, triggering a model update as soon as either method detects amplitude drift.

[0138] Compared to existing technologies that rely on a single indicator to determine model degradation, this invention provides a detection method sensitive to abrupt amplitude drift by dividing the relative deviation sequence into two equal segments and comparing the difference between the means of the two segments. This complements the linear regression slope method, enabling the model update mechanism to respond to both gradual and abrupt amplitude drift modes simultaneously.

[0139] Application examples

[0140] In a water diversion tunnel project, an open-face TBM was used for tunneling. The TBM had a cutterhead diameter of 7830 mm, a rated torque of 8800 kNm, a cutterhead rotation speed range of 0 to 7.6 rpm, a maximum thrust of 27488 kN, and a maximum advance speed of 120 mm / min. The tunnel traversed a geological formation primarily composed of granite, with uniaxial compressive strength varying between 80 MPa and 150 MPa. The TBM's control computer used a programmable logic controller (PLC) as its core processor, and the data acquisition module's sampling period was set to 1 second.

[0141] During tunneling, the industrial control computer collects the three-phase current signal of the cutterhead drive motor through a current transformer installed on the cutterhead drive motor, with a sampling frequency set to 2 kHz; it collects the propulsion pressure signal through a pressure sensor installed at the oil inlet of the propulsion hydraulic cylinder, also with a sampling frequency set to 2 kHz; it collects the cutterhead rotation speed through a rotary encoder installed on the cutterhead drive shaft; it obtains the propulsion force by combining the pressure sensor installed at the oil inlet of the propulsion hydraulic cylinder with the effective working area of ​​the hydraulic cylinder; it calculates the cutterhead torque by using the current and rotation speed of the cutterhead drive motor; and it collects the body attitude angle through an angle sensor installed on the TBM body, and uses the deviation angle of the body attitude angle relative to the design axis as the body attitude angle deviation value. At the end of each sampling period, the five parameters—cutterhead rotation speed, cutterhead torque, propulsion speed, propulsion force, and body attitude angle deviation value—are packaged in a fixed data frame format and stored in the historical tunneling dataset in chronological order.

[0142] After the TBM enters the stable tunneling phase, the industrial control computer reads 30 data frames from the historical tunneling dataset prior to the current time window, extracts five feature values ​​within that window: propulsion speed, mean propulsion force, propulsion force variance, fuselage attitude angle deviation, and rate of change of fuselage attitude angle deviation, and inputs these values ​​into a pre-trained Long Short-Term Memory (LSTM) network model. After forward computation, the model outputs predicted propulsion force and torque values, and accordingly determines the propulsion force and torque adjustments. In sections with relatively homogeneous lithology, such as granite with a uniaxial compressive strength of approximately 100 MPa, the model's output propulsion force adjustment is typically between -200 kN and +300 kN, and the torque adjustment is typically between -50 kNm and +80 kNm.

[0143] Meanwhile, within each acquisition cycle, a linear variable differential transformer displacement sensor installed on the electro-hydraulic proportional valve acquires the valve core displacement signal in real time, and a wire-type displacement sensor installed on the piston rod of the propulsion hydraulic cylinder acquires the piston displacement signal in real time. Both acquire the signal synchronously with a sampling period of 1 second and store it in a storage queue. The length of this storage queue is pre-determined through a step response test: a step control signal is sent to the electro-hydraulic proportional valve, and the time required from the signal transmission to the piston displacement reaching 90% of its steady-state value is recorded at a frequency of not less than 100 Hz. This time was measured to be approximately 4.5 seconds. Dividing this by the sampling period of 1 second and rounding up, the length of the storage queue is set to 5 sampling periods, meaning the delay time of the propulsion hydraulic system is approximately 5 sampling periods. When the storage queue is full, all the data in the queue are arranged in reverse order from latest to earliest acquisition time to form an inversion response sequence. Its length corresponds to the delay time and is not affected by the number of time steps within the judgment period.

[0144] A judgment cycle is triggered after each preset number of cutterhead rotations. In homogeneous strata, the preset number of rotations is set to 3. Taking a typical judgment cycle as an example, the cutterhead speed is 4 revolutions per minute, one sampling cycle is 1 second, and the time required for 3 rotations is approximately 45 seconds. At each time step within this judgment cycle, the industrial control computer reads the actual value of the hydraulic cylinder piston displacement at that time step from the historical tunneling dataset, and simultaneously reads the historical piston displacement value at a historical moment that is one time delay ahead from the historical tunneling dataset. The two values ​​are subtracted and the absolute value is taken to obtain the displacement deviation value for that time step. At the end of the judgment cycle, the displacement deviation values ​​of all time steps are summed and divided by the total number of time steps to obtain the average displacement deviation. If the preset maximum deviation is set to 2 mm and the measured average displacement deviation is 0.6 mm, then the ratio of the average displacement deviation to the preset maximum deviation is 0.3, the first weight is equal to 0.7, and the second weight is equal to 0.3. The valve core displacement signal at the current time step is multiplied by 0.7 and then added to the valve core displacement signal at the corresponding time in the inversion response sequence, which is multiplied by 0.3. The result of the addition is used as the valve core displacement compensation signal output.

[0145] The update ratio of the storage queue is determined based on the same average displacement deviation. The historical average displacement deviation of this TBM under normal operating conditions is 0.5 mm. The first deviation threshold is set at 0.4 mm, and the second deviation threshold is set at 1.0 mm. The measured average displacement deviation in this judgment period is 0.6 mm, which falls between the first deviation threshold of 0.4 mm and the second deviation threshold of 1.0 mm; therefore, the contents of the storage queue remain unchanged. If the measured average displacement deviation reaches 1.2 mm in a sudden change zone in hard rock, which is greater than the second deviation threshold of 1.0 mm, the storage queue is updated at the second update ratio of 80%. That is, 80% of the data stored in the storage queue is replaced. The amount of data replaced is the product of the storage queue length (5) and the 80% ratio, which is 4 sampling periods. When updating at the first update ratio of 20%, the amount of data replaced is 1 sampling period.

[0146] The propulsion adjustment output from the Long Short-Term Memory (LSTM) network model is converted from digital to analog to generate a propulsion control signal. This signal is the target displacement command for the electro-hydraulic proportional valve spool, measured in millimeters. This target displacement command is added to the spool displacement compensation signal output from step four, and the sum is sent as the final drive signal to the power amplifier of the electro-hydraulic proportional valve, driving the spool to the target position. The torque adjustment output from the LTM network model is also converted from digital to analog to generate a torque control signal, which is sent to the frequency converter of the cutterhead drive system. The frequency converter adjusts its output frequency and voltage based on the received signal, thereby regulating the output torque of the cutterhead drive motor.

[0147] During tunneling, the Long Short-Term Memory (LSTM) network model is updated online in the following manner: An update cycle is defined as 10 rotations of the cutterhead. At the end of each update cycle, the predicted thrust and torque sequences output by the model at each time step within that cycle, along with the corresponding actual thrust and torque sequences, are collected. The sign agreement rate between the predicted thrust difference sequence and the actual thrust difference sequence is calculated and denoted as the thrust trend agreement rate; the sign agreement rate between the predicted torque difference sequence and the actual torque difference sequence is calculated and denoted as the torque trend agreement rate. Both the first and second agreement rate thresholds are set to 0.80. During a crossing of a fault fracture zone, the direction of change in the model's predicted values ​​begins to lag behind the direction of change in the actual values, and the thrust trend agreement rate drops to 0.65, below the first agreement rate threshold of 0.80. The industrial control computer adds all data frames from that update cycle as new training samples to the training set, performs incremental training on the model, and updates the model parameters. After the update, the model's prediction trend agreement rate for the fault fracture zone conditions recovers to above 0.85.

[0148] Assuming the trend consistency rate meets the standard, the industrial control computer further detects amplitude drift. Within one update cycle, the thrust trend consistency rate is 0.88, which is greater than the first consistency rate threshold of 0.80. The relative deviation between the predicted thrust value and the actual thrust value at each time step within this cycle is calculated to obtain the thrust relative deviation sequence. The linear regression slope of this sequence is calculated. If the slope is positive and greater than the preset positive slope threshold of 0.02, a positive amplitude drift is determined to have occurred, triggering incremental training. In addition, the relative deviation sequence is divided into a first segment and a second segment in chronological order, and the mean of the first segment and the mean of the second segment are calculated separately. If the difference between the mean of the second segment and the mean of the first segment is greater than the first difference threshold of 0.10, a positive abrupt amplitude drift is determined to have occurred, also triggering incremental training. The two amplitude drift detection methods complement each other, which is beneficial for timely model updates in both gradual and abrupt amplitude drift scenarios.

[0149] The above steps are continuously cycled during the tunneling process, realizing intelligent control of TBM construction.

[0150] The number of devices and processing scale described herein are for the purpose of simplifying the description of the invention. Applications, modifications, and variations of the invention will be readily apparent to those skilled in the art.

[0151] Although embodiments of the present invention have been disclosed above, they are not limited to the applications listed in the specification and embodiments. They can be applied to various fields suitable for the present invention. For those skilled in the art, other modifications can be easily made. Therefore, without departing from the general concept defined by the claims and their equivalents, the present invention is not limited to the specific details and embodiments shown and described herein.

Claims

1. A TBM construction intelligent control method, characterized in that, Includes the following steps: Step 1: Collect cutterhead rotation speed, cutterhead torque, propulsion speed, propulsion force, and machine attitude angle deviation values ​​at a fixed sampling period, and package them into data frames and store them in the historical tunneling dataset; Step 2: Input the propulsion speed, mean propulsion force, variance of propulsion force, fuselage attitude angle deviation and rate of change of the historical tunneling dataset from 10 to 50 data frames before the current time step into the long short-term memory network model to obtain the propulsion force adjustment and torque adjustment. Step 3: Collect valve core displacement signal and piston displacement signal and store them in historical tunneling dataset and storage queue. The length of the dataset is equal to the number of sampling cycles from receiving the command from the valve core to the piston producing displacement. After the dataset is full, arrange them in reverse order to form the inversion response sequence. Step 4: Using a preset number of rotations of the cutterhead as the judgment period, calculate the average displacement deviation between the actual piston displacement value within this judgment period and the historical piston displacement value in the historical tunneling dataset one delay time earlier. The ratio of this average displacement deviation to the preset maximum deviation is used as the second weight. When the value of the second weight is greater than 1, it is taken as 1. 1 minus the second weight is used as the first weight. The current valve core displacement signal and the valve core displacement signal at the corresponding time in the inversion response sequence are weighted and summed according to the two weights to obtain the valve core displacement compensation signal. When the average displacement deviation is less than the first deviation threshold, the storage queue is updated with the first update ratio; when it is greater than the second deviation threshold, the storage queue is updated with the second update ratio. Other cases remain unchanged. The second deviation threshold is greater than the first deviation threshold. Step 5: Convert the thrust adjustment amount into a valve core target displacement command and superimpose it with the compensation signal before sending it to the electro-hydraulic proportional valve; convert the torque adjustment amount into a torque control signal and send it to the frequency converter of the cutter head drive system.

2. The intelligent control method for TBM construction as described in claim 1, characterized in that, In step one, the current signal of the cutter head drive motor and the pressure signal of the propulsion hydraulic system are also collected simultaneously. Between step one and step two, the following steps are also included, wherein one collection cycle is defined as one revolution of the cutter head: Envelope detection is performed on the current signal acquired in step one to obtain the current amplitude envelope signal. Fourier transform is performed on the current amplitude envelope signal to obtain the envelope spectrum. The fundamental frequency amplitude and the second harmonic amplitude of the cutter head rotation frequency are extracted from the envelope spectrum. The ratio of the second harmonic amplitude to the fundamental frequency amplitude is calculated and recorded as the second harmonic ratio. For the pressure signal of the propulsion hydraulic system collected in step one, calculate the number of rising edges of the pressure signal per unit time, and record the number of rising edges as the pressure impact frequency. The second harmonic ratio is multiplied by the pressure impact frequency to obtain the geological characteristic fusion value; When the geological feature fusion value increases monotonically over three consecutive acquisition cycles, it is determined that the rock at the working face has hardened. The mean thrust extracted in step two is multiplied by a first correction coefficient greater than 1 and then used as the corrected mean thrust input into the long short-term memory network model. When the geological feature fusion value decreases monotonically over three consecutive acquisition cycles, the lithology of the tunnel face is determined to be softened. The mean thrust extracted in step two is multiplied by a second correction coefficient less than 1 and then used as the corrected mean thrust input into the long short-term memory network model.

3. The intelligent control method for TBM construction as described in claim 2, characterized in that, The current signal of the cutterhead drive motor and the pressure signal of the propulsion hydraulic system are collected. The current signal collected in step one is subjected to envelope detection to obtain the current amplitude envelope signal. The specific method for obtaining the envelope spectrum by performing Fourier transform on the current amplitude envelope signal is as follows: One acquisition cycle is defined as one revolution of the cutter head, and current and pressure signals are continuously acquired within each acquisition cycle. The current signal of the current acquisition cycle and the previous four consecutive acquisition cycles, a total of five acquisition cycles, are stored in the current sliding window; the pressure signal of the current acquisition cycle and the previous two consecutive acquisition cycles, a total of three acquisition cycles, are stored in the pressure sliding window. At the end of each acquisition cycle, the envelope spectrum is obtained by performing envelope detection and Fourier transform on all current signals within the current sliding window. The number of rising edges is counted for all pressure signals within the pressure sliding window. This number is divided by the total duration of the pressure sliding window to obtain the pressure impact frequency of the pressure sliding window. The geological feature fusion value for that acquisition period is obtained by multiplying the second harmonic ratio of the current sliding window with the pressure impact frequency of the pressure sliding window.

4. The intelligent control method for TBM construction as described in claim 3, characterized in that, At the end of each acquisition cycle, the envelope spectrum is obtained by performing Fourier transform on all current signals within the current sliding window and counting the number of rising edges on all pressure signals within the pressure sliding window. This number is then divided by the total duration of the pressure sliding window to obtain the pressure impact frequency of that pressure sliding window. At the end of each acquisition cycle, the intersection of the current sliding window and the pressure sliding window on the time axis is defined as the common time window, and the length of the common time window is equal to the length of the pressure sliding window. After envelope detection of all current signals within the current sliding window, Fourier transform is performed to obtain the envelope spectrum, and the second harmonic ratio is extracted from the envelope spectrum. The number of rising edges of the pressure signal within the common time window is counted, and this number is divided by the total duration of the common time window to obtain the pressure impact frequency of the common time window. The geological feature fusion value for that acquisition period is obtained by multiplying the second harmonic ratio by the pressure impact frequency of the common time window.

5. The intelligent control method for TBM construction as described in claim 1, characterized in that, The specific method for obtaining the thrust adjustment and torque adjustment from the long short-term memory network model in step two is as follows: The Long Short-Term Memory network model outputs a thrust prediction and a torque prediction at each time step; The instantaneous thrust error is obtained by subtracting the predicted thrust value output at the current time step from the actual thrust value collected at the current time step, and the instantaneous torque error is obtained by subtracting the predicted torque value output at the current time step from the actual torque value collected at the current time step. The instantaneous thrust errors of P consecutive time steps, including the current time step, are used to form a thrust error sequence. The variance of the thrust error sequence is calculated and denoted as the thrust error variance. The instantaneous torque errors of P consecutive time steps, including the current time step, are used to form a torque error sequence. The variance of the torque error sequence is calculated and denoted as the torque error variance. When the thrust error variance is greater than the first variance threshold, the thrust prediction value of the current time step is replaced with the arithmetic mean or weighted moving average of the actual thrust values ​​of the previous preset number of time steps and then output as the thrust adjustment amount. The preset number is less than P. When the torque error variance is greater than the second variance threshold, the torque prediction value of the current time step is replaced with the arithmetic mean or weighted moving average of the actual torque values ​​of the previous preset number of time steps, and then used as the torque adjustment amount output.

6. The intelligent control method for TBM construction as described in claim 5, characterized in that, When the thrust error variance exceeds the first variance threshold, the specific method for replacing the thrust prediction value at the current time step with the arithmetic mean or weighted moving average of the actual thrust values ​​at the previous preset number of time steps as the thrust adjustment output is as follows: Calculate the variance of the propulsion error sequence for P consecutive time steps, including the current time step, and denote it as the current variance; For each of the P consecutive time steps including the current time step, calculate the difference between the predicted thrust value of the time step and the actual thrust value offset by k time steps before the time step, where k ranges from 1 to K and K is the preset maximum delay step number, and obtain K offset error sequences. Calculate the variance of each offset error sequence and denote it as the offset variance. The current variance is compared with each offset variance one by one. When there is an offset variance smaller than the current variance, the k value corresponding to the smallest offset variance is selected as the estimated delay step; when each offset variance is greater than or equal to the current variance, the estimated delay step is 0. When the estimated delay step is 0 or 1, the predicted thrust value of the current time step is replaced with the actual thrust value of the previous time step and then output as the thrust adjustment amount. When the estimated delay step number is greater than 1 and less than or equal to the preset first delay threshold, the predicted thrust value of the current time step is replaced with the arithmetic mean of the actual thrust values ​​of the estimated delay step number of time steps before the current time step, and then the thrust adjustment amount is output. When the estimated delay step count is greater than the first delay threshold, the predicted thrust value at the current time step is replaced with the weighted moving average of the actual thrust values ​​at a preset fixed number of time steps prior to the current time step, and this is used as the thrust adjustment amount output.

7. The intelligent control method for TBM construction as described in claim 6, characterized in that, The specific methods for calculating the variance of the thrust error sequence over P consecutive time steps, including the current time step, and for calculating the variance of each offset error sequence are as follows: Perform linear regression on the thrust error sequence to obtain the linear regression slope of the thrust error sequence; A trend term sequence is generated based on the slope and intercept of the linear regression. The trend value at the corresponding time in the trend term sequence is subtracted from each value in the thrust error sequence to obtain the detrended thrust error sequence. Calculate the variance of the detrended driving force error sequence and use it as the current variance; After performing linear detrending on each offset error sequence, the variance of each offset error sequence is calculated and used as the offset variance.

8. The intelligent control method for TBM construction as described in claim 1, characterized in that, In step two, the long short-term memory network model is updated online in the following manner: An update cycle is defined as N revolutions of the cutterhead, where N is an integer greater than 1. At the end of each update cycle, the propulsion prediction sequence and torque prediction sequence output by the long short-term memory network model at each time step within that update cycle are collected, along with the actual propulsion sequence and actual torque sequence at the corresponding time. Calculate the difference between the predicted values ​​of two adjacent time steps in the thrust prediction value sequence to obtain the thrust prediction difference sequence. Calculate the difference between the actual values ​​of two adjacent time steps in the actual thrust value sequence to obtain the actual thrust difference sequence. Calculate the sign consistency rate between the thrust prediction difference sequence and the actual thrust difference sequence, denoted as the thrust trend consistency rate. Calculate the difference between the predicted values ​​of two adjacent time steps in the torque prediction value sequence to obtain the torque prediction difference sequence. Calculate the difference between the actual values ​​of two adjacent time steps in the actual torque value sequence to obtain the actual torque difference sequence. Calculate the sign consistency rate between the torque prediction difference sequence and the actual torque difference sequence, denoted as the torque trend consistency rate. When the thrust trend consistency rate is lower than the first consistency rate threshold or the torque trend consistency rate is lower than the second consistency rate threshold, all data frames within the update period are added as new training samples to the training set of the long short-term memory network model. The model parameters are then updated after incremental training of the long short-term memory network model.

9. The intelligent control method for TBM construction as described in claim 8, characterized in that, It also includes the following steps: When the thrust trend consistency rate is greater than or equal to the first consistency rate threshold, the relative deviation between the thrust prediction value and the actual thrust value at each time step in the update cycle is further calculated. The relative deviation is the thrust prediction value minus the actual thrust value and divided by the non-zero actual thrust value to obtain the thrust relative deviation sequence. Calculate the slope of the linear regression of the relative deviation sequence of propulsion force; When the slope of the linear regression is positive and greater than the preset positive slope threshold, or when the slope of the linear regression is negative and the absolute value is greater than the absolute value of the preset negative slope threshold, all data frames in the update period are added as new training samples to the training set of the long short-term memory network model, and the model parameters are updated after incremental training of the long short-term memory network model.

10. The intelligent control method for TBM construction as described in claim 9, characterized in that, Mutant amplitude drift was also detected in the following manner: Calculate the relative deviation between the predicted thrust value and the actual thrust value at each time step within the update cycle, arrange them in chronological order to obtain the thrust relative deviation sequence, and divide the sequence into a first segment sequence and a second segment sequence in chronological order. When the sequence length is odd, the middle element is included in the first segment sequence. Calculate the arithmetic mean of the first part of the sequence as the first part mean, and calculate the arithmetic mean of the second part of the sequence as the second part mean. Calculate the difference between the mean of the second segment and the mean of the first segment; When the difference is greater than the first difference threshold, or when the difference is less than the second difference threshold, all data frames in the update period are added as new training samples to the training set of the long short-term memory network model. The model parameters are updated after incremental training of the long short-term memory network model. The second difference threshold is negative.