Wave-flow coupling scouring simulation test method in super-gravity field

Abstract

The invention discloses a wave-current coupling scouring simulation test method in a high gravity field, and relates to the crossing field of geotechnical engineering and ocean engineering. According to the method, multi-source observation signals are synchronously collected in a supergravity environment, and a phase synchronous observation sequence with a wave phase as an index is constructed. A low-dimensional filtering-inversion integrated method based on wave phase harmonic unwrapping is utilized, and a bed surface shear stress field is inverted with high precision in combination with an adaptive basis function and scouring starting physical constraint. And carrying out shear stress controllability analysis by calculating the sensitivity of the shear stress to the wave current phase difference and constructing a nonlinear controllable domain. Based on this, a double-layer control strategy is adopted, an outer layer plans a target phase trajectory, an inner layer is combined with a bed surface evolution state to switch a linear or nonlinear compensation control structure, and a control instruction for driving a wave generation and current generation mechanism is generated. According to the invention, the problems that the dynamic change of the wave flow phase cannot be simulated and the real-time feedback of the shear stress field is lacked in the prior art are solved, and the precise closed-loop simulation of the wave flow coupling scouring is realized.

Description

Technical Field

[0001] This invention belongs to the interdisciplinary technical field of geotechnical engineering and marine engineering, and in particular to the simulation and testing method of wave-current coupling scour in hypergravity fields. Background Technology

[0002] Offshore engineering structures such as port projects, cross-sea bridges, and subsea pipelines face severe wave and current scouring challenges during long-term service, and their foundation stability is directly related to national energy and transportation security. In-depth understanding of the wave-current coupling scouring mechanism in complex marine environments and calculation of the evolution of bed shear stress are of significant strategic importance for engineering safety assessment and life prediction.

[0003] Currently, research on localized scour mainly relies on conventional gravity field (1g) model experiments or early centrifuge model experiments. In existing centrifuge simulation techniques, wave generation and current generation usually operate independently, or can only achieve a simple linear superposition of waves and water flow. For experiments involving wave-current coupling, existing equipment mostly uses mechanical locking to maintain a fixed wave-current phase relationship, such as simulating only the same or opposite directions, and bed surface monitoring mainly relies on sparsely arranged point sensors to obtain single-point time history data of flow velocity or pressure.

[0004] However, existing technologies suffer from deep technical bottlenecks in simulating real marine environments and achieving precise control, including distortion in the simulation of wave-current dynamic coupling effects and a lack of bed surface response feedback control. Specifically: the fixed-phase loading mode severs the spatiotemporal connection between the wave-current phase difference and the dynamic evolution of tides, making it impossible to reproduce the nonlinear modulation process of shear stress by wave-current interaction; traditional discrete-point observation and time-domain filtering methods are difficult to isolate the non-steady-state interference of bed surface morphology (such as sand waves) evolution on the flow field, resulting in insufficient accuracy of shear stress field inversion and an inability to construct an effective mapping relationship between phase difference and shear stress; due to the lack of real-time assessment of shear stress controllability, the control system struggles to adaptively adjust the phase adjustment strategy when faced with bed surface response lag or nonlinear saturation, causing the experimental process to deviate from the preset scour evolution trajectory. Summary of the Invention

[0005] The purpose of this invention is to provide a method for simulating wave-current coupling scour in a hypergravity field, in order to solve one of the problems mentioned above in the existing technology.

[0006] Technical solution: A method for simulating and testing wave-current coupling scour in a hypergravity field, comprising:

[0007] Read the monitoring signals from multiple sensors and the wave-generating and current-generating motion state signals under the hypergravity test environment, analyze the real-time wave phase and water flow phase, map the monitoring signals to the wave phase domain, and generate a phase synchronization observation sequence;

[0008] Based on the phase-synchronous observation sequence, the fluid dynamic characteristics of the bed surface are inverted and analyzed to obtain the estimated data of shear stress field reflecting the stress state of the bed surface;

[0009] Based on the deviation between the estimated shear stress field data and the preset target shear stress index, and combined with the sensitivity characteristics of the shear stress field to changes in wave-current phase difference, wave-current phase control commands are generated for real-time adjustment of the phase difference between the wave-generating mechanism and the current-generating mechanism.

[0010] Beneficial effects: This invention solves the problems of existing technologies being unable to simulate dynamic changes in wave-current phase and lacking real-time feedback of shear stress field, and realizes accurate closed-loop simulation of wave-current coupled scouring. Attached Figure Description

[0011] Figure 1 This is a flowchart illustrating the steps of a wave-current coupling scour simulation test method in a hypergravity field according to an embodiment of this application.

[0012] Figure 2 This is a flowchart illustrating the steps involved in generating a phase-synchronized observation sequence in an embodiment of this application.

[0013] Figure 3 This is a flowchart illustrating the steps involved in calculating and evaluating the scour response in an embodiment of this application.

[0014] Figure 4 This is a flowchart illustrating the steps for generating the shear stress-phase difference correlation sequence and calculating the shear stress sensitivity index in this embodiment of the application.

[0015] Figure 5 This is a three-dimensional structural diagram of the centrifugal simulation system for the scouring process in a hypergravity field in this embodiment of the application.

[0016] Figure 6 This is a top view of the centrifugal simulation system for the scouring process in a hypergravity field, as described in this application.

[0017] Figure 7 For along Figure 6 Schematic diagram of the cross-sectional structure along line AA in the middle.

[0018] Figure 8 For along Figure 6 Schematic diagram of the cross-sectional structure of the middle BB line.

[0019] In the diagram: 1. Wave height phase sensor; 2. Displacement sensor; 3. Flow meter; 4. Wave damping net; 5. Annular water tank; 6. Water tank support structure; 7. Servo cylinder; 8. Transparent plexiglass; 9. High-precision camera; 10. Water flow phase sensor; 11. Phase output terminal; 12. Phase control transmission mechanism; 13. Phase input terminal; 14. Miniature pore water pressure sensor; 15. Laser profilometer; 16. Hydraulic model; 17. Honeycomb flow stabilizer; 18. Shear stress sensor; 19. Pulsating pressure sensor; 20. Sand; 21. Soil sample box; 22. Propeller. Detailed Implementation

[0020] Combination Figures 1 to 8 This describes the implementation details of this application.

[0021] Example 1: Construct a centrifugal simulation system that can accurately realize wave-current coupled loading, measure bed surface shear stress in real time, and capture local scouring process in a hypergravity field, solving the problem that existing equipment cannot realize dynamic phase control of wave-current and full-field measurement of shear stress.

[0022] The wave-current coupling scour simulation system mainly consists of three core subsystems: a dynamically phase-adjustable wave-current coupling system, a multi-source data acquisition system, and a bed shear stress inversion and control system. The wave-current coupling system generates waves and water flow within the model chamber at the end of the centrifuge arm and dynamically adjusts the phase between them. The multi-source data acquisition system synchronously acquires various physical quantities during the experiment. The bed shear stress inversion and control system processes the data and generates feedback control commands.

[0023] Specifically, the dynamically phase-adjustable wave-current coupling system includes a wave-generating device, a current-generating device, and a phase control transmission mechanism. The wave-generating device employs a wedge-shaped wave generator. Driven by a servo cylinder, the wedge reciprocates vertically, generating waves by displacing water. A wave-damping net is installed at the end. The wedge's motion period is adjustable from 0 to 10 seconds, capable of generating regular and irregular waves with heights of 0 to 10 meters. The current-generating device uses a propeller, capable of generating constant or pulsating flow at speeds of 0 to 3 m / s. The phase control transmission mechanism includes a harmonic reducer and a gear linkage device. The phase input is connected to the servo drive shaft of the wave-generating wedge, and the phase output is connected to the motor shaft of the current-generating propeller. Through a combination of mechanical linkage and electronic servo control, the phase difference between wave generation and water flow pulsation is adjusted, achieving dynamic phase changes from 0 to 180 degrees.

[0024] In terms of the multi-source data acquisition system, this embodiment is equipped with a hydrodynamic sensor group, a scour sensor group, and a bed surface observation sensor group. The hydrodynamic sensor group includes a wave height phase sensor located upstream of the test section to measure wave surface elevation changes; a current meter located near the bottom to measure the water flow velocity vector; and a pulsating pressure sensor embedded in the model surface to acquire the dynamic pressure distribution under wave-current action. The scour sensor group includes a non-contact laser displacement sensor for real-time monitoring of bed surface scour depth; and a micropore water pressure sensor embedded inside the sand sample box to monitor seepage and excess pore pressure response within the sand. The bed surface observation sensor group includes a particle image velocimetry (PIV) system for acquiring near-bottom flow field structure and a digital image correlation system for monitoring bed surface micro-deformation. In the test environment, a hydraulic model simulating structures such as piles and piers was installed to study local scour effects. A honeycomb flow stabilizer was installed to ensure the quality of the incoming flow and stabilize the flow field. All sensors are triggered by a high-precision synchronous clock, and the acquired data has a unified time reference.

[0025] The bed shear stress inversion and control system, serving as the core of data processing and decision-making, is deployed in the host computer workstation. This system receives data from multiple sensor sources, inverts the bed shear stress field using a built-in algorithm, and calculates control commands based on preset test targets, such as the target shear stress time history or target scour depth. Through a servo controller, it adjusts the motion parameters of the wave-generating wedge and the rotational speed parameters of the current-generating propeller, achieving closed-loop control of the wave-current coupling conditions.

[0026] In this embodiment, to adapt to the hypergravity environment, the model box and all water-contact components are made of high-strength, corrosion-resistant materials, such as stainless steel or aluminum alloy. Sensor signal transmission uses shielded cables resistant to electromagnetic interference, fixed along the inner wall of the model box to reduce interference with the flow field. During centrifuge operation, the gravitational acceleration is set to N times the gravitational acceleration, and the model's geometric dimensions are scaled down by 1 / N to meet the similarity criteria for centrifuge model experiments.

[0027] Example 2 describes how to accurately extract the phase information of waves and water flow from the original sensor signals in a supergravity wave-flow coupling experiment, and construct a synchronous observation sequence indexed by the wave phase.

[0028] Step 201: Read the monitoring signals from the multi-source sensors and the wave-generating and current-generating motion state signals under the hypergravity test environment, and analyze the real-time wave phase and water flow phase.

[0029] Specifically, during the experiment, the data acquisition system received signals in real time from various monitoring devices such as wave height phase sensors, flow meters, and pressure sensors, as well as status feedback signals from the wave generator and flow generator servo controllers. To obtain an accurate phase reference, this embodiment employs specific signal processing methods to extract the wave phase and water flow phase separately.

[0030] For real-time wave phase analysis, the wave height phase sensor signal located upstream of the test section is extracted from the multi-source sensor monitoring signals. Since the wave surface exhibits periodic fluctuations over time, the peak detection algorithm can be used to identify the wave peak time. Let t... wave_peak Let Tw be the most recently detected wave peak time, and Tw be the current wave period, which can be calculated in real time by the time difference between adjacent wave peaks or given by the wave generation command. Then, the real-time wave phase θ at the current time t is... w It can be calculated using linear interpolation, with the formula: θ w = 360° * (t - t) wave_peak ) / Tw. This phase value cycles between 0 and 360°, defining the peak as 0° or 360°. To improve noise immunity, the wave height signal can be low-pass filtered before peak detection to remove high-frequency glitches.

[0031] The analysis of real-time water flow phase primarily addresses situations involving pulsating flow. This involves extracting the propeller rotation speed pulse signal or encoder feedback signal from the wave-generating motion state signal. Since there is a positive correlation between water flow velocity and propeller rotation speed: v = k * n, the pulsation period Tf of the flow velocity is typically determined by the propeller rotation speed adjustment period. This is achieved by identifying the peak time t of the rotation speed or flow velocity signal. flow_peak And the pulsation period Tf can be used to calculate the real-time water flow phase θ at the current time t. f The formula is: θ f = 360° * (t - t) flow_peak ) / Tf. If it is a steady flow, the flow phase can be considered constant or not included in the phase difference calculation.

[0032] Based on this, a high-precision synchronous clock signal is used to align the real-time wave phase and the real-time water flow phase with timestamps. Since different sensors and control systems may have different sampling rates, the system uses a unified master clock, such as a PTP protocol or GPS clock, as a reference, and all data points are timestamped to the microsecond level. By aligning the timestamps, the wave-current phase difference Δ at the same moment can be accurately calculated. θ = θ w - θ f .

[0033] Step 202: Map the monitoring signal to the wave phase domain to generate a phase-synchronous observation sequence.

[0034] This step is used to convert time series data, which originally had time t as the independent variable, into data with wave phase θ as the independent variable. w The phase domain data is for the independent variable.

[0035] Based on the real-time wave phase, the time-continuous multi-source sensor monitoring signal is divided into several complete wave cycle units.

[0036] Specifically, the system detects θ w The moment when the wave transitions from 360° back to 0° or crosses 0° is marked as the start of a new wave cycle. The data segment between two adjacent zero-phase moments constitutes a wave cycle unit. Each unit contains time-series segments of all sensors within that cycle, such as wave height, flow velocity, pressure, and shear stress.

[0037] For wave periodic units, the instantaneous signal values ​​at preset discrete phase nodes are extracted to form a periodic segmented observation set.

[0038] In this embodiment, the preset discrete phase nodes can be uniformly distributed, for example, taking one point every 10°, such as 0°, 10°, ..., 350°, for a total of 36 nodes. For each node, the data point corresponding to the phase position is found in the wave cycle unit. If the sampling time of the original data does not fall exactly on that type of phase node, linear interpolation or spline interpolation methods are used to calculate the instantaneous signal value at that phase node. For the k-th wave cycle, the observation vector D is obtained. k (θ), where θ takes all preset phase nodes.

[0039] By statistically superimposing the periodic segmented observation set at the same discrete phase node, a phase-synchronized observation sequence with wave phase as the independent variable index is generated.

[0040] Specifically, since a single measurement signal may contain random turbulent fluctuations and measurement noise, this embodiment employs a phase averaging method to extract the deterministic wave period component. The data D from N consecutive wave periods is then processed. k The phase-synchronized observation sequence D is obtained by superimposing and averaging (θ) (k=1...N). avg (θ) = (1 / N) * Σ(D k (θ)). In addition, statistics such as standard deviation can be calculated to characterize the stability of turbulence intensity or periodic variations.

[0041] In some alternative implementations, to address minor fluctuations in wave cycles, a dynamic time warping algorithm can be used to align data from different cycles before averaging. For non-stationary processes, such as those with varying wave heights or periods, a sliding window averaging method can be used, which only superimposes the most recent M cycles to preserve the time-varying characteristics of the process.

[0042] It should also be noted that the generated phase-synchronized observation sequence not only includes the mean values ​​of the original physical quantities, but also constructs a bidirectional index mapping structure between the time domain and the phase domain. This index table records which original time points the phase domain data points originated from and their corresponding wave cycle IDs. This bidirectional indexing mechanism allows for rapid backtracking to the original time series for verification if anomalies are detected in the phase domain data during subsequent analysis, or for accurately mapping the phase domain inversion results back to the time axis for real-time control.

[0043] Example 3 describes the harmonic inversion and multi-source fusion process based on physical constraints.

[0044] Based on the phase-synchronous observation sequence, the fluid dynamic characteristics of the bed surface are inverted and analyzed to obtain the estimated data of shear stress field reflecting the stress state of the bed surface.

[0045] This embodiment details a shear stress field inversion method based on wave phase harmonic expansion. Utilizing a phase-synchronous observation sequence, it constructs a harmonic model incorporating physical constraints to achieve high-precision, low-dimensional reconstruction of the bed surface shear stress field. Specifically:

[0046] Step 301: A pre-defined shear stress harmonic model is used to describe the change of shear stress at the bed surface spatial nodes with wave phase. The model uses the harmonic coefficient vector and phase parameter vector as the state variables to be solved.

[0047] In this step, based on the periodic characteristics of wave action, it is assumed that the shear stress τ(x,t) at any point x on the bed surface can be represented as the superposition of the fundamental wave and its higher harmonics. Since the data has been mapped to the phase domain, the independent variable here is the wave phase θ. The shear stress harmonic model is specifically expressed as: τ(x, θ) = a0(x) + Σ m=1 M [a m (x) * cos(m * θ +φ m (x))] where a0(x) represents the average shear stress, mainly contributed by steady-state water flow, M is the cutoff order, for example, 3 or 5, a m (x) is the amplitude coefficient of the m-th harmonic, φ m (x) represents the phase parameter of the m-th harmonic. The state variable to be solved is the parameter vector S(x) = [a0, a1, ..., a...] at each spatial node. M, φ1, ..., φ M This modeling approach compresses the high-dimensional shear stress field, which originally changes rapidly over time, into a low-dimensional parameter field that changes slowly over space, thereby reducing the dimensionality and computational cost of the inversion problem.

[0048] Step 302: Analyze the characteristics of bed elevation change in the multi-source sensor monitoring signals, determine the bed evolution stage marker corresponding to the current monitoring area, the marker includes the flat bed stage or the sand wave development stage, and select the harmonic basis function family accordingly.

[0049] To improve the model's adaptability, this embodiment introduces an adaptive basis function selection mechanism. The bed surface elevation data z obtained using a laser profilometer or multibeam echo sounder is used. bed (x, t) is used to calculate its spatial spectrum or spatial gradient. If the root mean square height of the bed elevation is less than the preset threshold and there is no obvious main peak in the spatial spectrum, it is determined to be the flat bed stage; if there are spatial periodic fluctuations and the waveform shows obvious asymmetry, with the upstream side being gentle and the downstream side being steep, it is determined to be the sand wave development stage.

[0050] When the bed evolution stage is designated as the flat bed stage, a symmetric family of harmonic basis functions consisting of sine and cosine functions is selected. At this point, the standard Fourier series expansion described above is sufficient to describe the symmetric periodic variation of the shear stress.

[0051] When the bed evolution stage is marked as the sand wave development stage, an asymmetric harmonic basis function family containing asymmetric waveform characteristic functions is selected. For example, second-order skew terms can be introduced into the basis functions, or non-sinusoidal basis functions such as Jacobi elliptic functions can be used to better capture the shear stress abrupt changes and asymmetries caused by the flow separation on the back surface of the sand wave.

[0052] The selected family of harmonic basis functions is used to construct the functional expression of the shear stress harmonic model.

[0053] Step 303: Establish an optimization solution system with the objective of minimizing the difference between the output of the shear stress harmonic model and the phase-synchronous observation sequence, and introduce a physical constraint term for scour.

[0054] To solve for the above model parameters, the following optimization objective function J is constructed: J = J data + λ * J smooth +γ *J physics Among them, J data This is the data fidelity term, representing the difference between the shear stress calculated by the model and the measured or extrapolated shear stress in the phase-synchronous observation sequence. It is usually expressed in the least squares form: Σ ||G(S) - D obs || 2G is the observation operator, which maps the shear stress parameters to the observation space, such as the velocity gradient measured by PIV. D obs This is the actual observation sequence. smooth For spatial smoothing regularization, used to suppress non-physical spatial oscillations, for example: Σ ||Laplacian(S)|| 2 λ is the regularization parameter.

[0055] Specifically, this embodiment introduces J physics As a physical constraint term for scouring initiation, a critical shear stress threshold τ is introduced to describe the initiation conditions of the sediment used in the experiment. cr This threshold can be determined using a Shields curve based on the sediment particle size.

[0056] Based on the shear stress amplitude calculated from the critical shear stress threshold and the shear stress harmonic model, a scour-induced physical constraint term is constructed. The scour-induced physical constraint term is configured to accumulate numerical values ​​only for spatial distribution regions where the shear stress amplitude exceeds the critical shear stress threshold.

[0057] Specifically, it can take the form of a penalty function: J physics = ∫ [ max(0, |τ(x, θ)| - τ Max_physical ) 2 [dx] or constraints related to starting conditions, the shear stress in the region where scour has not occurred should not exceed τ. cr In this embodiment, the constraint term is preferably used to penalize regions where extremely high shear stress is physically impossible, such as non-scour hotspots, thereby improving the physical plausibility of the inversion results. The scour-induced physical constraint term is added as a weighted penalty component to the objective function of the optimization solution system, and the weighting coefficient γ can be adjusted according to the signal-to-noise ratio.

[0058] Step 304: Calculate the optimal harmonic coefficient vector and phase parameter vector using the optimization solution system, and reconstruct the shear stress field estimation data with unified spatiotemporal resolution based on these vectors.

[0059] Nonlinear optimization algorithms, such as the Levenberg-Marquardt algorithm or sequential quadratic programming, are used to solve the objective function J to obtain the optimal state parameter vector S. opt (x). After obtaining the optimal parameters, substituting them into the harmonic model formula, the shear stress value τ(x, θ) at any phase θ and any position x can be reconstructed. Since there is a clear mapping relationship between θ and time t, shear stress field estimation data τ(x, t) with uniform spatiotemporal resolution can be further generated.

[0060] In some preferred embodiments, to further improve accuracy, especially for local high-gradient regions such as the edges of scour pits, a multi-source data fusion strategy is employed. Traditional inversion methods, such as those based on extended Kalman filtering or simple interpolation, use the shear stress field obtained as the basis, i.e., the low-frequency component. High-frequency detail features, such as higher-order harmonic components, are extracted using the parameters obtained from the harmonic inversion described above, and the two are then superimposed.

[0061] Specifically, the operation can be as follows: in the identified high shear stress gradient regions, the weight of higher harmonics is increased; while in flat regions, the fundamental wave or traditional inversion results are mainly relied upon. The final output shear stress field retains the overall macroscopic trend while clearly depicting the fine structure of local scour hotspots.

[0062] Example 4: Construction of shear stress controllability analysis and nonlinear controllable domain.

[0063] Based on the deviation between the estimated shear stress field data and the preset target shear stress index, and combined with the sensitivity characteristics of the shear stress field to changes in wave-current phase difference, wave-current phase control commands are generated for real-time adjustment of the phase difference between the wave-generating mechanism and the current-generating mechanism.

[0064] Due to the evolution of the bed surface morphology, such as the generation and movement of sand waves and the nonlinear hydrodynamic effects under hypergravity, the influence of the wave-current phase difference on the bed surface shear stress is not constant, but rather time-varying and state-dependent. Therefore, it is necessary to assess in real time whether adjusting the phase difference is effective under the current operating conditions and to what range it should be adjusted.

[0065] Step 401: Generate the shear stress-phase difference correlation sequence and calculate the shear stress sensitivity index.

[0066] The key shear stress index sequence, used as the control object, is extracted from the shear stress field estimation data. These key shear stress indices are scalars or vectors characterizing the current scour state, and may include: the spatially average effective shear stress within the monitoring area, the maximum peak shear stress, or the local shear stress at specific scour hotspot locations, such as around a structure. A wave-current phase difference sequence is obtained by using real-time wave phases and real-time water flow phases and calculating the difference between them. The shear stress index sequence and the wave-current phase difference sequence are aligned on the time axis to form a correlation sequence.

[0067] Based on this, a sliding time window technique is used to fit the local relationship between the key shear stress index sequence and the wave-current phase difference sequence, and the shear stress sensitivity index reflecting the shear stress response intensity caused by the unit phase difference change is calculated. A sliding window containing several wave cycles is set, and the rate of change of the shear stress index with respect to the wave-current phase difference is fitted within the window. This embodiment proposes a multi-timescale sensitivity analysis mechanism:

[0068] Based on a preset timescale of one or a few wave cycles (e.g., 1-3 cycles) in the multi-source sensor monitoring signals, linear regression or difference methods are used to fit the rapid response relationship between the key shear stress index sequence and the wave-current phase difference sequence, thus obtaining the short-term instantaneous sensitivity κ. short It is used to quickly adjust the phase difference to achieve immediate control of shear stress.

[0069] Based on the pre-defined timescale of bed morphology evolution in the phase-synchronous observation sequence, such as dozens of cycles or the characteristic time of sand wave migration, the long-term cumulative trend of the key shear stress index sequence as a function of the wave-current phase difference sequence is fitted to obtain the mid-term evolution trend parameter κ. Medium This is used to assess the long-term impact of phase difference adjustment on the direction of bed evolution.

[0070] Step 402: Construct a controllable domain for the nonlinear wave-current phase difference.

[0071] To overcome the limitations of traditional linear control in the strongly nonlinear region, this embodiment introduces the concept of a nonlinear wave-current phase difference controllable domain. This controllable domain Ω(t) is defined as the set of all possible wave-current phase difference values ​​that can keep the predicted shear stress index within a preset target range, i.e., a safe range, at the current time t.

[0072] Based on short-term instantaneous sensitivity and medium-term evolution trend parameters, the possible wave-current phase difference range is scanned and predicted at the current time section to determine the set of values ​​that can make the predicted shear stress value fall within the target range, thus forming a nonlinear wave-current phase difference controllable domain.

[0073] The specific construction process is as follows: A shear stress harmonic model is used as a fast predictor. At the current time section, keeping the bed surface morphology and fluid parameters constant, a virtual scan is performed on the possible range of wave-current phase difference Δθ from 0 to 360° or from 0 to 2π radians. For the virtual phase difference values ​​within the scan range, they are substituted into the harmonic model to calculate the corresponding predicted shear stress index. The predicted results are then compared with the upper limit τ of the target shear stress index. Max and lower limit τ Min Compare the predicted values. If the predicted value falls within [τ]... Min , τ Max If the virtual phase difference is within the interval, then the virtual phase difference is included in the controllable domain Ω(t).

[0074] By combining the numerical range of the target shear stress index with the shear stress sensitivity index, the effective adjustment range of the wave-current phase difference that can keep the shear stress within a preset safe range is determined, generating shear stress controllability analysis data. The final generated shear stress controllability analysis data includes not only the aforementioned sensitivity index κ and controllable domain Ω(t), but also a controllability level determination. For example, when Ω(t) is an empty set, it is determined to be uncontrollable, and the system may need to alarm or adjust other parameters such as wave height and flow rate; when the Ω(t) range is wide and the sensitivity κ is moderate, it is determined to be highly controllable, suitable for fine-tuning.

[0075] In some alternative implementations, the sensitivity index can be calculated using the recursive least squares method to estimate the sensitivity coefficient online, and a forgetting factor can be introduced to adapt to the time-varying characteristics of the system. When constructing the controllable domain, physical constraints of the actuator, such as the upper limit of the phase adjustment speed, can also be considered. Phase differences that cannot be reached within the next control cycle are removed from the controllable domain to obtain a dynamically reachable controllable domain.

[0076] Example 5 describes the specific implementation method of wave-current phase closed-loop control using controllability analysis results. This example adopts a two-layer control structure and introduces a controller structure switching mechanism based on the physical state of the bed surface to solve the control problems caused by shear stress response lag and nonlinear enhancement during sand wave development.

[0077] Step 501, outer control loop: calculate the target wave flow phase difference trajectory.

[0078] The outer control loop is primarily responsible for path planning. Within the outer control loop, the deviation e between the current estimated shear stress field and the target shear stress index is calculated. τ The controllability analysis data of shear stress and the controllable domain of nonlinear wave-current phase difference are used as constraints. Based on the deviation e... τ The sign and short-term instantaneous sensitivity κ short The sign of the variable determines the direction of phase difference adjustment. For example, if the current shear stress is too high and the sensitivity is positive, the phase difference should be reduced.

[0079] Next, the target wave current phase difference trajectory Δθ is calculated within the next control cycle. target The calculation logic is: Δθ target = Δθ current + α * (e τ / κ short ), where α is the outer layer adjustment gain. To ensure safety, Δθ is calculated. targetThe control must be confined within the controllable domain Ω(t) of the nonlinear wave-current phase difference. If the calculated value exceeds the boundary of Ω(t), the boundary value is taken as the target, or the optimal operating point within Ω(t) is searched based on the mid-term evolution trend parameters. Throughout this process, the control target remains within the physically safe region.

[0080] Step 502, Inner control loop: Generation of instructions based on bed surface state switching.

[0081] In the inner control loop, the actual phase difference is obtained by calculating the difference between the real-time wave phase and the real-time water flow phase. Based on the tracking error of the target wave-current phase difference trajectory according to the actual phase difference, wave-current phase control commands are generated to drive the wave-generating mechanism and the current-generating mechanism to perform phase adjustment actions. The inner control loop is responsible for high-frequency tracking of the target trajectory given by the outer layer. The system analyzes the bed surface morphology estimation data and determines the current bed surface evolution stage marker in real time by calculating the spatial spectrum characteristics of the bed surface elevation, such as wavelength, wave height, or topographic gradient. The inner control loop is equipped with a controller structure switching mechanism based on the physical state of the bed surface. It analyzes the spatial spectrum characteristics and topographic gradient in the bed surface morphology estimation data to determine the current bed surface evolution stage marker, which includes the flat bed stage and the sand wave development stage.

[0082] When the bed evolution stage is marked as the flat bed stage, the response of shear stress to the phase difference is approximately linear and hysteresis-free. At this time, the linear phase tracking control structure is activated, and the wave-current phase control command is calculated using a linear feedback algorithm to achieve rapid tracking of the target wave-current phase difference trajectory. This structure uses an incremental PID control algorithm. Specifically, the difference between the target wave-current phase difference trajectory and the actual phase difference at the current sampling time k is calculated to obtain the phase tracking error e(k). The feedback adjustment amount Δu(k) is calculated using the proportional coefficient Kp, integral coefficient Ki, and derivative coefficient Kd. The specific formula is: Δu(k) = Kp * [e(k) - e(k-1)] + Ki * e(k) + Kd * [e(k) - 2e(k-1) + e(k-2)]. The first derivative of the target trajectory is calculated and multiplied by the feedforward coefficient kf to generate the feedforward compensation amount u. ff The final wave-current phase control command is the sum of the feedback adjustment and the feedforward compensation.

[0083] When the bed surface evolution stage is marked as the sand wave development stage, the shear stress response caused by the shear stress due to the vortex shedding on the back surface of the sand wave exhibits a nonlinear hysteresis effect, making it difficult for a simple linear PID controller to maintain control quality. At this point, the system activates a nonlinear compensation control structure, superimposing a nonlinear compensation term on top of the linear feedback algorithm to offset the asymmetric shear stress response hysteresis caused by the sand wave terrain, generating wave-current phase control commands. Based on the aforementioned linear PID algorithm, a nonlinear compensation term u is superimposed... nlThis compensation term is constructed based on the asymmetric coefficients in the shear stress harmonic model, such as the second harmonic amplitude, and is used to pre-compensate for shear stress waveform distortion caused by sand wave topography. During this stage, the system dynamically adjusts the PID parameters according to the sensitivity index; for example, it appropriately increases Kp to enhance tracking capability when sensitivity decreases.

[0084] In some preferred embodiments, the inner control loop also includes a dead zone. When the phase tracking error e(k) is less than a preset dead zone threshold, such as a phase difference of 1°, no adjustment is performed to avoid mechanical wear caused by frequent fine-tuning of the actuator. Before the control commands are output to the wave-generating and current-generating mechanisms, they undergo amplitude limiting and rate limiting processing, ensuring that the commands remain within the physical safety operating range of the equipment.

[0085] Example 6 describes how, after completing a single-round wave-current coupling scour test, multi-source data is used to calculate and evaluate the test results, and the test parameters for the next round are optimized accordingly, forming a large closed loop of test-evaluation-optimization.

[0086] Step 601, Calculate and evaluate the scour response.

[0087] Based on wave height data, flow velocity data, and pulsating pressure data from multi-source sensor monitoring signals, the shear stress under the action of waves alone, the shear stress under the action of water flow alone, and the wave-flow coupled interactive shear stress are calculated respectively. The effective shear stress is obtained by superimposing the shear stress under the action of waves alone, the shear stress under the action of water flow alone, and the wave-flow coupled interactive shear stress.

[0088] After executing one round of wave-current phase control commands, the system collects shear stress field estimation data, bed scour depth monitoring data, and pulsating pressure data for the entire process. It calculates the ratio of the effective shear stress to the sum of the shear stress caused by the wave alone and the shear stress caused by the water flow alone, generating a coupling effect coefficient used to calculate the degree of wave-current synergistic enhancement or attenuation.

[0089] To calculate the wave-current coupling effect, this embodiment uses the coupling effect coefficient η for evaluation. The calculation formula is: η = τ couple / (τ wave + τ current ), where τ couple τ is the measured effective shear stress. wave and τ current These are the shear stress components calculated separately based on wave and current theories, respectively. If η > 1, it indicates that wave-current coupling promotes scouring; if η < 1, it indicates that scouring is inhibited.

[0090] Establish a computational correlation model for the scouring rate. Based on the physical model in the disclosure document, the scouring rate v... erode It can be represented as: v erode = k * τcouple α * ΔP β Where k, α, and β are undetermined coefficients, and ΔP is the pulsating pressure. The scour depth variation curve over time, measured in this round of experiments, is used, i.e., v... erode The measured values ​​and the inverted shear stress and pressure data were fitted with the coefficients k, α, and β using multivariate nonlinear regression. These coefficients reflect the unique scouring behavior of sandy materials under hypergravity conditions.

[0091] Step 602: Parameter iteration and next round of target generation.

[0092] Based on the above evaluation results, scour response evaluation data is generated. If the scour effect of this round of tests does not achieve the expected failure mode, such as insufficient scour pit depth or excessive shear stress fluctuations during the control process, the system will correct the load control target data for the next round of tests.

[0093] Specific correction strategies include: deriving the ideal shear stress trajectory required to reach the target scour depth based on the fitted scour rate model, and updating the target shear stress index data. Utilizing controllable domain information, the change pattern of the target wave-current phase difference is redesigned. For example, if a phase difference interval, such as around 90 degrees, is found to have the largest coupling coefficient η and high control sensitivity in the previous test, the duration of this phase difference interval will be increased in the next test to accelerate the scour process.

[0094] The updated test load control target data, including the new target wave sequence, flow sequence, and phase difference planning, is used as input for the next round of testing, and the online control process is restarted. Through multiple rounds of iterative approximation, accurate simulation of the local scour mechanism of complex marine structures and accurate determination of the ultimate bearing capacity can be achieved.

[0095] In summary, this invention presents a complete simulation technology for supergravity wave-current coupling scour, from low-level hardware perception to high-level algorithm decision-making. It solves the problem of unclear visualization by using phase synchronization and harmonic inversion, and solves the problem of uncontrollability by using controllability analysis and dual-layer switching control, thus realizing a true closed-loop simulation of wave-current coupling.

[0096] Example 7 describes the rapid estimation and multimodal fusion inversion of instantaneous shear stress.

[0097] Although harmonic inversion methods have advantages in handling periodic wave and current characteristics, in the short timescale of real-time control, or in the initial loading stage where harmonic characteristics have not yet stabilized, instantaneous estimation based on physical approximation and state fusion based on traditional filtering are still the key means to ensure the robustness of the system.

[0098] Step 701: Generate initial calculation data for instantaneous shear stress.

[0099] To meet the low latency requirements of real-time control, typically millisecond-level latency is required. This step employs an explicit formula based on the boundary layer physical approximation to directly calculate the instantaneous shear stress. The system extracts the vertical component of flow velocity and the horizontal component of water pressure in the near-bed region from multi-source sensor monitoring signals. Specifically, flow velocity data at locations very close to the bed surface, such as 0.5 mm to 2 mm in height, is acquired using a PIV system or a high-frequency point-type velocimeter. Combined with this measurement point height, the near-bed flow velocity gradient du is calculated. dz The pressure gradient dp is calculated using an array of pulsating pressure sensors embedded in the bed surface along the flow direction. dx .

[0100] Based on the boundary layer momentum equation approximation in fluid dynamics, the initial calculation data τ0 of the instantaneous shear stress reflecting the transient forces acting on the bed surface is obtained. The calculation formula is: τ0 = μ * du dz + k p * dp dx Where μ is the dynamic viscosity coefficient of the fluid, and for sediment-laden water flows, this coefficient needs to be corrected according to the sediment concentration; k p This is the pressure gradient correction coefficient, used to characterize the nonlinear contribution of the pressure gradient to the bed surface shear stress. This coefficient was determined through previous clean water calibration tests. This calculation process does not involve complex iterative solutions, therefore it is extremely fast and can serve as a rapid feedback signal in the real-time control loop.

[0101] Step 702: Construct a multimodal fusion filtering system.

[0102] To overcome the problems of high noise and insufficient spatial coverage of single sensors, this embodiment constructs a multimodal fusion filtering module. This module establishes a unified state vector and observation equations. The state vector s(t) is defined to include three core physical quantities: the bed shear stress field τ(x, t), the near-bed velocity profile characteristics u(z, t), and the bed elevation z. bed (x, t). That is, s(t) = {τ, u, z} bed}

[0103] Constructing nonlinear state transition equations s k +1 = f(s k ) + w k The function f(s) k This describes the evolution of the system state over time, utilizing the sediment initiation and transport equations to describe the coupling relationship between shear stress and bed elevation changes. k Characterize the uncertainty and noise of the model. Construct the observation equation d. k = H k * s k + vk Wherein, the observation vector d k It includes real-time readings from all sensors, such as PIV flow rate, pressure, laser displacement, and the observation matrix H. k This describes how state variables are mapped to the observation space, v k To observe noise.

[0104] Step 703: Perform extended Kalman filtering and regularized inversion.

[0105] The extended Kalman filter algorithm is used to recursively solve the above state-space model. At each time step, a time update is performed based on the state estimate and state transition equation from the previous time step, i.e., the prediction step, to obtain the prior state estimate; the current observation vector d is then used... k The measurement is updated using the observation equation, i.e., the correction step, and the posterior state estimate τ is obtained by weighting the residuals with Kalman gain. hat (x, t) and the corresponding error covariance matrix P(x, t).

[0106] In some preferred embodiments, to obtain a higher spatial resolution shear stress field, the system periodically, for example, runs a regularized high-resolution inversion module every 10 time steps. This module uses the posterior estimate of the extended Kalman filter output as initial values ​​to construct an optimization problem that includes data fidelity terms and spatial smoothing regularization terms: min || G * τ - d|| 2 + λ * || L * τ || 2 Where G is the forward operator, L is the Laplacian operator or gradient operator, and λ is the regularization parameter, for example, 0.01. This inversion process uses the alternating direction multiplier method or the conjugate gradient method to solve the problem, and outputs a full-field shear stress distribution with a spatial resolution of centimeters.

[0107] Example 8: A centrifugal simulation system for scouring processes in a hypergravity field, consisting of three parts: a dynamically phase-adjustable wave-current coupling system, a multi-source data acquisition system, and a bed surface shear stress inversion system.

[0108] The dynamically phase-adjustable wave-current coupling system includes: a propeller, a wedge block, a phase sensor, a transmission device, and a servo control device.

[0109] The propeller can generate water flow of 0-3 m / s, and the wedge block can generate regular and irregular waves with a period of 0-10 s and a wave height of 0-10 m. The phase sensors include wave height phase sensors and water flow phase sensors.

[0110] The wave phase sensor, employing a miniature level gauge, is non-contact to avoid interfering with the flow field. It is installed 10cm upstream of the test section in the model box, with the probe facing the water surface. The water flow phase sensor is installed at the shaft end of the flow-generating propeller motor. It obtains the water flow phase by detecting rotational speed pulses. The flow velocity is positively correlated with the propeller speed, and the rotational speed phase is the water flow phase.

[0111] Wave phase extraction method: Using the wave height-time curve output by the level gauge, the peak time t is identified by a peak detection algorithm. 波峰 trough moment t 波谷 Define wave phase θ w =360°×(tt 波峰 ) / T w T w For the wave period, θ w The range is 0-360°.

[0112] Water flow phase extraction method: Calculate the water flow velocity u=k based on the propeller rotation speed n (r / min). u ×n,k u Define the flow phase θ as the velocity-rotation speed calibration coefficient. f =360°×(tt 流速峰值 ) / T f T f The period of water flow pulsation is T; when there is no pulsation, T... f Take 10 times the wave period to ensure phase synchronization.

[0113] The transmission device includes a harmonic reducer and a gear linkage mechanism. The input end is connected to the wave-generating wedge block servo device to obtain the displacement signal of the wedge block's up and down movement. The output end is connected to the flow-generating propeller speed shaft through the gear linkage mechanism to realize wave and flow linkage.

[0114] The wave-generating wedge is driven by a servo hydraulic cylinder, and the current-generating propeller is driven by a permanent magnet synchronous motor.

[0115] The specific steps of the closed-loop control algorithm for real-time phase adjustment, which uses the host computer software to preset arbitrary phase difference variation patterns, are as follows:

[0116] Phase difference calculation: Real-time acquisition of θ w and θ f Calculate Δθ 实际 =θ w -θ f If the value is negative, add 360° to convert it to the range of 0-360°;

[0117] Error handling: Calculate the deviation e = Δθ 目标 -Δθ 实际 Range-extended PID control is adopted:

[0118] Δu(k)=K p [e(k)-e(k-1)]+K i e(k)+K d [e(k)-2e(k-1)+e(k-2)], where K p =0.5, K i =0.01, K d =0.1;

[0119] Feedforward compensation: based on the first derivative dt / Δθ of the preset phase difference curve. 实际 / dt, output adjustment amount in advance, u 前馈 =k f ×dt / Δθ 实际 / dt, where k f =0.2;

[0120] Execution output: Total adjustment u(k) = Δu(k) + u 前馈 This is converted into control signals for wave generation and current generation;

[0121] Wave generation: the velocity of the wedge block, u w =k w ×u(k), where k w =0.1mm / s per;

[0122] Current generation: propeller speed n f =n f0 +k f ×u(k), where k f =5r / min per°, n f0 The reference speed;

[0123] The signal flow path is: phase sensor → data acquisition card → host computer → servo controller (executes PID algorithm) → wave / current generation mechanism → phase sensor (feedback closed loop). By using a synchronized clock, the timestamp error between sensor acquisition, algorithm calculation, and actuator action is ensured to be ≤1ms, achieving wave-current synchronization.

[0124] Dynamically phase-tunable wave-current coupling mechanism:

[0125] Wave-current motion linkage is achieved through a phase sensor and a transmission device: the wave phase sensor identifies the wave crest moment in real time, and the water flow phase sensor detects the peak flow velocity synchronously. The two are mechanically linked with the transmission device through a harmonic reducer to ensure the timing of the motion is matched.

[0126] The preset phase difference curve adjusts the time difference between the wedge block and the propeller through a PID algorithm, achieving stepless adjustment from 0 to 180° with a phase error of ≤ ±3°, accurately simulating the wave and current phase relationship of the real ocean.

[0127] The multi-source data acquisition system consists of a hydrodynamic sensor group, a local scour sensor group, and a bed surface observation sensor group.

[0128] The hydrodynamic sensor set includes:

[0129] Wave height phase sensor: One sensor is installed upstream (10cm from the hydraulic model) and one downstream (10cm from the hydraulic model) of the test section, 5cm above the water surface, measuring vertically downwards;

[0130] Flow meter: Two flow meters are installed on the test section, one upstream and one downstream, 5cm from the bed surface (wave-current coupling core area), with the probes facing the direction of water flow;

[0131] Shear stress sensors: Three sensors are directly attached to the soil sample box bed surface, 6cm in front of the hydraulic model, and evenly distributed.

[0132] Pulsating pressure sensors: One sensor is embedded in the flow-facing surface of the hydraulic model (one at the bottom, one in the middle, and one at the top). The sensor probe is flush with the surface of the model to avoid interfering with the flow field.

[0133] Effective shear stress τ couple =τ wave +τ current +τ interact , where τ wave The shear stress (Pa) caused by the wave acting alone is calculated from wave height phase sensor data, τ wave =0.5×ρ×g×H² / (2πL), where ρ is the density of water (1000kg / m³), g is the gravitational acceleration, H is the wave height, and L is the wavelength;

[0134] τ current The shear stress (Pa) acting solely on the water flow is calculated from the flow meter data, τ. current =f×ρ×v² / 8, where f is the bed surface friction coefficient, calibrated by the baseline data of the shear stress sensor; v is the average flow velocity of the water.

[0135] τ interact The shear stress (Pa) of wave-current coupling is derived from data from a pulsating pressure sensor, τ. interact =k1×ΔP, where k1 is the interaction coefficient, ranging from 0.05 to 0.1; ΔP is the pulsating pressure.

[0136] The flushing sensor group includes:

[0137] Displacement sensor: mounted directly above the viewing window on the top of the model box, aligned with the center of the soil sample box test section, with the laser vertically irradiating the soil surface to achieve non-contact measurement of penetration depth;

[0138] Pore ​​water pressure sensor: buried inside the soil (one at a depth of 5cm and one at a depth of 10cm), 10cm away from the edge of the hydraulic model foundation to avoid the sensor damaging the soil structure, and outputs the bottom pore pressure.

[0139] Calculation of coupling effect coefficient η: η = τ couple / (τ wave +τ current );

[0140] η>1: Wave-current synergistic enhancement (coupling promotes scouring);

[0141] η=1: No synergistic effect (wave and current act independently);

[0142] η<1: Wave-current synergistic weakening (coupling suppresses scouring).

[0143] The bed surface observation sensor group includes:

[0144] Near-bottom velocity field acquisition module (PIV): Acquires instantaneous velocity u(z,x,t) at a distance of 0.3-2mm from the bed surface through a thin-film laser and a high-speed camera, and outputs a near-bed velocity profile u(z)=Ψu / Ψz, where Ψ is the partial derivative and z is the height from the bed surface;

[0145] Digital Imaging (DIC) System: A flexible film is laid on the bed surface, and the micro-deformation ε(x,t) of the bed surface is acquired by a high-speed imaging system;

[0146] The bed surface observation sensor group is integrated into a high-precision camera, and all sensors in the system achieve the same source of data acquisition through a high-precision synchronous clock.

[0147] The system is externally equipped with an annular water tank, a water tank support structure, and transparent acrylic glass. The annular water tank is installed within the main frame of the system via the water tank support structure, and its side walls are made of transparent acrylic glass, which facilitates non-contact measurement and recording of the wave-current coupling scouring process in the hypergravity field by a high-precision camera.

[0148] Bed surface morphology measurement module: Uses a laser profilometer to scan and obtain the bed surface elevation z. bed (x,t);

[0149] Wire treatment: Electromagnetic interference shielded wires are used, with a wire diameter ≤1mm, to reduce interference with the flow field. The wires are laid along the inner wall of the model box, and the outlet is sealed with sealant.

[0150] Synchronous acquisition is achieved using a synchronous clock with a set acquisition frequency of 100Hz. The synchronous clock sends a unified trigger signal, and all sensors start acquiring data simultaneously, avoiding time differences caused by asynchronous acquisition.

[0151] Before the experiment, each sensor was calibrated in a 1g environment, the calibration coefficients were recorded, and the data were written into the acquisition system.

[0152] Start the centrifuge to the required hypergravity acceleration for the experiment, and collect baseline data from the sensors under undulating loading conditions to ensure data accuracy.

[0153] scouring rate v erode Model associated with coupling parameters: v erode =k×τ couple α ×ΔP β ;

[0154] Among them, v erode The instantaneous scouring velocity (cm / s) is calculated from displacement sensor data: v erode =(h t -h {t-1} ) / Δt,h t For the depth of the charge at time t, h {t-1} The depth of the charge at time t-1 is Δt=1ms. k, α, β: scaling correction coefficients for the hypergravity field, calibrated experimentally.

[0155] Specifically, the core of wave-current coupling effect quantification is to transform physical formulas into programming logic and then into real-time calculations. The algorithm is embedded in the servo control system and linked with the acquisition module to achieve real-time acquisition of one set of data, calculation of one result, and output of one quantification index.

[0156] The coupling effect quantification algorithm of this invention is based on multi-field synchronous data acquisition, and the specific steps are as follows:

[0157] Step 1, Data Preprocessing:

[0158] Acquire synchronously the wave height H, average flow velocity v, pulsating pressure ΔP, and depth h. t Data such as wave height H and flow velocity v, with timestamp errors ≤ 1ms, and piercing depth h. t Record continuously at a sampling interval Δt = 1ms;

[0159] Step 2, Calculation of individual shear stress:

[0160] a. Wave-induced shear stress τ wave The wavelength L is obtained by fitting data from the wave height phase sensor, and then calculated using the formula τ. wave =0.5×ρ×g×H² / (2πL) is used for calculation, where ρ is the density of water (1000 kg / m³), g is the gravitational acceleration, and under hypergravity, it is replaced by n×g, where n is the centrifugal coefficient;

[0161] b. Shear stress τ under the sole action of water flow current According to formula τ current= f × ρ × v² / 8, where f is the coefficient of friction of the bed surface, calibrated using baseline data from the shear stress sensor, with a calibration error ≤ 5%;

[0162] Step 3, Calculation of Coupled Interactive Shear Stress:

[0163] Based on the ΔP (unit kPa) collected by the pulsating pressure sensor, according to the formula τ interact =k1×ΔP×1000 (ΔP converted to Pa), where k1 is the interaction coefficient of the hypergravity field, obtained by fitting experimental data, with a fitting range of 0.05-0.1 and a fitting error ≤3%;

[0164] Step 4: Calculation of effective shear stress and coupling effect coefficient:

[0165] a. Effective shear stress τ couple =τ wave +τ current +τ interact ;

[0166] b. Coupling effect coefficient η=τ couple / (τ wave +τ current ), and truncate outliers (η∈[0.8,1.5]);

[0167] Step 5: Establish a quantitative correlation model between effective shear stress, pulsating pressure, and erosion rate:

[0168] According to formula v erode =(h t -h {t-1} The measured erosion rate was calculated using Δt, and then a correlation model was established using the hypergravity scaling correction factors k, α, and β. erode =k×τ couple α ×ΔP β ;

[0169] Step 6, Real-time Correction:

[0170] Compare the measured flushing rate v every 1 minute. erode Compared with the predicted value v erode_pred If the error exceeds 10%, the coefficients of k, α, and β are refitted to ensure calculation accuracy.

[0171] The bed shear stress inversion system is the core module of data processing. It inverts the spatiotemporal field τ(x,t) from multi-source observations d(t) and provides uncertainty estimates to meet real-time control requirements. The core principle of the inversion algorithm includes three levels: a fast instantaneous shear stress estimation module, a multi-modal fusion filtering module, and a regularized high-resolution inversion.

[0172] The bed shear stress inversion system is connected in series with the multi-source data acquisition system to obtain synchronous observation data.

[0173] a. Instantaneous shear stress rapid estimation module:

[0174] Based on a combined approximation of velocity and pressure gradients, it features low latency (1-3ms) and is used for real-time control: τ0(x,t)=μ*Ψu / Ψz| z=0 +k p *Ψp / Ψx, where Ψ is the partial derivative and μ is the dynamic viscosity.

[0175] b. Multimodal fusion filtering module:

[0176] Establish the state vector: s(t)=[τ(x,t),u(z,t),z bed [(x,t)], where τ(x,t) is the shear stress field, u(z,t) is the near-bottom velocity profile, and z bed (x,t) represents the bed surface elevation;

[0177] Construct the state equation: s k+1 =f(s k )+w k , where f(s) k ) is the state prediction function, w k Characterizing model uncertainty;

[0178] Construct the observation equation: d k =H k *s k +v k , where d k H represents the observation vector. k Represents the multimodal observation matrix, v k Indicates observation noise;

[0179] The extended Kalman filter (EKF) is used for real-time fusion, and the output is the inverted shear stress τ̂(x,t) and the covariance matrix P(x,t).

[0180] c. Regularized high-resolution inversion module:

[0181] Solve using regularization with bed shape evolution as a constraint: min||G*τ-d|| 2 +λ||L*τ|| 2 , where G is the positive operator matrix, d is the observation data vector, λ=0.01 is the regularization parameter, L is the regularization operator, outputs the full-field shear stress distribution τ(x,t) with a spatial resolution of 1cm×1cm, and the calculation delay is ≤1s.

[0182] The specific implementation steps are as follows:

[0183] a. System alignment and time synchronization: Start PPS and confirm that the timestamps of all acquisition terminals are consistent (error <1ms).

[0184] b. Sensor zero-point and gain calibration: Measure and record the zero point and noise spectrum under static conditions.

[0185] c. Thin film calibration: Apply a known force / displacement and establish a calibration curve G for the thin film displacement → apparent τ.

[0186] d. Pore pressure → τ calibration: Apply a known τ to the circulating water tank using a reference wall shear stress sensor and record the pore pressure response.

[0187] e. Velocity gradient path calibration: The correction coefficient μ*Ψu / Ψz is obtained by synchronously measuring with PIV and comparing with the reference shear stress, where Ψ is the partial derivative.

[0188] f. Matrix construction: Summarize all observed data into an H matrix.

[0189] g. Start the instantaneous shear stress fast estimation module: output τ0 for each frame and send it to the controller.

[0190] h. Parallel execution of the multimodal fusion filtering module: Extended Kalman filter (EKF) fuses observations from 20-200Hz to update the inverted shear stress τ̂ and covariance P.

[0191] i. Periodically execute the regularized high-resolution inversion module: perform high-precision inversion on the data of the most recent time window using the alternating direction multiplier method, and save the results for analysis.

[0192] j. In each control loop, τ̂ is sent back to the host computer, and the phase module makes corresponding closed-loop adjustments.

[0193] To address the issue of simulation distortion in the dynamic coupling effect of waves and current, a dynamically adjustable phase hardware system combined with a two-layer wave-current control strategy was adopted. The outer control loop plans a phase difference trajectory that conforms to real tidal patterns, while the nonlinear compensation mechanism of the inner control loop overcomes the hysteresis effect of sand wave topography on the flow field. This achieves continuous and precise dynamic adjustment of the wave-current phase difference within the range of 0-180°, thus reproducing the real physical process of wave-current coupling.

[0194] To address the issues of missing feedback control for bed response and insufficient accuracy in shear stress field inversion, a shear stress field inversion method based on wave phase harmonic expansion was adopted. This method abandons traditional discrete point measurements and utilizes phase-synchronous observation sequences and adaptive harmonic basis functions to reconstruct a continuous shear stress field from sparse observations, effectively eliminating interference from bed morphology evolution. Simultaneously, by introducing shear stress controllability analysis and nonlinear controllable domain construction, a mapping between phase adjustment and shear stress response was established, enabling the control system to possess the ability to perceive controllability. This allows the control strategy to adaptively adjust according to the real-time bed state, ensuring the scouring process strictly follows a preset trajectory.

[0195] In summary, the instantaneous estimation and extended Kalman filtering fusion method provided by this invention provides high-quality initial guesses for harmonic inversion, accelerating the convergence of harmonic parameters. On the other hand, it directly serves as a feedback source for the control system in transient processes where the harmonic model has not yet been established or the wave flow phase is unclear, ensuring the continuity and safety of control throughout the entire time period.

Claims

1. A method for simulating wave-current coupling scour in a hypergravity field, characterized in that, include: Read the monitoring signals from multiple sensors and the wave-generating and current-generating motion state signals under the hypergravity test environment, analyze the real-time wave phase and water flow phase, map the monitoring signals to the wave phase domain, and generate a phase synchronization observation sequence; Based on the phase-synchronous observation sequence, the fluid dynamic characteristics of the bed surface are inverted and analyzed to obtain the estimated data of shear stress field reflecting the stress state of the bed surface; Based on the deviation between the estimated shear stress field data and the preset target shear stress index, and combined with the sensitivity characteristics of the shear stress field to changes in wave-current phase difference, wave-current phase control commands are generated for real-time adjustment of the phase difference between the wave-generating mechanism and the current-generating mechanism.

2. The method according to claim 1, characterized in that, Mapping the monitoring signals to the wave phase domain to generate a phase-synchronized observation sequence includes: Based on the real-time wave phase, the time-continuous multi-source sensor monitoring signal is divided into several complete wave cycle units; For wave periodic units, the instantaneous signal values ​​at preset discrete phase nodes are extracted to form a periodic segmented observation set; By statistically superimposing the periodic segmented observation set at the same discrete phase node, a phase-synchronized observation sequence with wave phase as the independent variable index is generated.

3. The method according to claim 2, characterized in that, Based on the phase-synchronous observation sequence, the hydrodynamic characteristics of the bed surface were inverted and analyzed, and the estimated shear stress field data reflecting the stress state of the bed surface were obtained, including: A pre-defined shear stress harmonic model is used to describe the variation of shear stress at the bed surface nodes with wave phase. The model uses the harmonic coefficient vector and phase parameter vector as the state variables to be solved. An optimization solution system is established with the objective of minimizing the difference between the output of the shear stress harmonic model and the phase-synchronized observation sequence; The optimal harmonic coefficient vector and phase parameter vector are calculated using an optimization solution system, and the spatiotemporal resolution-uniform shear stress field estimation data are reconstructed accordingly.

4. The method according to claim 3, characterized in that, The pre-defined shear stress harmonic model describing the variation of shear stress at the bed surface nodes with wave phase includes: The characteristics of bed elevation change in multi-source sensor monitoring signals were analyzed to determine the bed evolution stage markers corresponding to the current monitoring area, including the flat bed stage or the sand wave development stage. When the bed surface evolution stage is marked as the flat bed stage, a symmetric harmonic basis function family composed of sine and cosine functions is selected; When the bed surface evolution stage is marked as the sand wave development stage, an asymmetric harmonic basis function family containing asymmetric waveform characteristic functions is selected; The selected family of harmonic basis functions is used to construct the functional expression of the shear stress harmonic model.

5. The method according to claim 3, characterized in that, An optimization solution system with the objective of minimizing the difference between the output of the shear stress harmonic model and the phase-synchronized observation sequence includes: A critical shear stress threshold is introduced to describe the initiation conditions of the sediment used in the experiment; Based on the shear stress amplitude calculated from the critical shear stress threshold and the shear stress harmonic model, a scour-induced physical constraint term is constructed. The scour-induced physical constraint term is configured to accumulate values ​​only for spatial distribution regions where the shear stress amplitude exceeds the critical shear stress threshold. The physical constraints caused by scouring are added as weighted penalty components to the objective function of the optimization solution system.

6. The method according to claim 1, characterized in that, Inverting and analyzing the fluid dynamics characteristics of the bed surface also includes calculating the coupling effect coefficients, and the specific steps include: Based on wave height data, flow velocity data and pulsating pressure data from multi-source sensor monitoring signals, the shear stress under the action of waves alone, the shear stress under the action of water flow alone, and the interactive shear stress of wave-flow coupling are calculated respectively. The effective shear stress is obtained by superimposing the shear stress caused by wave action alone, the shear stress caused by water flow alone, and the wave-flow coupled interactive shear stress. Calculate the ratio of the effective shear stress to the sum of the shear stress caused by the wave alone and the shear stress caused by the water flow alone, and generate a coupling effect coefficient for calculating the degree of wave-flow synergistic enhancement or attenuation.

7. The method according to claim 1, characterized in that, By combining the sensitivity characteristics of the shear stress field to changes in wave-current phase difference, shear stress controllability analysis data are generated. The specific steps include: The key shear stress index sequence, which is the control object, is extracted from the shear stress field estimation data, and the wave-current phase difference sequence is obtained by calculating the difference between the real-time wave phase and the real-time water flow phase. The sliding time window technique was used to fit the local relationship between the key shear stress index sequence and the wave-current phase difference sequence, and the shear stress sensitivity index reflecting the shear stress response intensity caused by the unit phase difference change was calculated. By combining the numerical range of the target shear stress index with the shear stress sensitivity index, the effective adjustment range of the wave-current phase difference that can keep the shear stress within the preset safe range is determined, and shear stress controllability analysis data is generated.

8. The method according to claim 7, characterized in that, The calculated shear stress sensitivity indices, which reflect the shear stress response intensity caused by a unit phase difference change, include: Based on the time scale of a single wave cycle preset in the multi-source sensor monitoring signal, the rapid response relationship between the key shear stress index sequence and the wave-current phase difference sequence is fitted to obtain the short-term instantaneous sensitivity. Based on the time scale of the bed surface morphology evolution preset in the phase synchronous observation sequence, the long-term cumulative trend of the key shear stress index sequence with the change of the wave-current phase difference sequence is fitted to obtain the mid-term evolution trend parameters. Based on short-term instantaneous sensitivity and medium-term evolution trend parameters, the possible wave-current phase difference range is scanned and predicted at the current time section to determine the set of values ​​that can make the predicted shear stress value fall within the target range, thus forming a nonlinear wave-current phase difference controllable domain.

9. The method according to claim 8, characterized in that, The wave-current phase control command generated for real-time adjustment of the phase difference between the wave-generating mechanism and the current-generating mechanism adopts a two-layer control structure, including: In the outer control loop, based on the deviation between the current shear stress field estimation data and the target shear stress index, and with the shear stress controllability analysis data and the nonlinear wave-current phase difference controllable domain as constraints, the target wave-current phase difference trajectory in the next control cycle is calculated. In the inner control loop, the difference between the real-time wave phase and the real-time water flow phase is calculated to obtain the actual phase difference. Based on the tracking error of the target wave-flow phase difference trajectory according to the actual phase difference, wave-flow phase control commands are generated to drive the wave-generating mechanism and the current-generating mechanism to perform phase adjustment actions.

10. The method according to claim 9, characterized in that, The inner control loop is equipped with a controller structure switching mechanism based on the physical state of the bed surface. The specific steps include: The spatial spectral characteristics and morphological gradients in the bed surface morphology estimation data are analyzed to determine the current bed surface evolution stage markers, including the flat bed stage and the sand wave development stage. When the bed evolution stage is marked as the flat bed stage, the linear phase tracking control structure is activated, and the wave-current phase control command is calculated using a linear feedback algorithm to achieve rapid tracking of the target wave-current phase difference trajectory. When the bed surface evolution stage is marked as the sand wave development stage, the nonlinear compensation control structure is activated. On the basis of the linear feedback algorithm, the nonlinear compensation term is superimposed to offset the asymmetric shear stress response lag caused by the sand wave terrain and generate wave flow phase control command.

Images