Open channel flow verification system based on attractor stability and multi-pump cooperation method

Abstract

The invention relates to an open channel flow verification system based on attractor stability and a multi-pump cooperation method. The open channel flow verification system comprises the steps that a reservoir, a water pump set with at least two water pumps connected in parallel, an open channel rectification section provided with a flow guide grid or a rectification grid, a flow stabilization section, a measurement section and a water outlet section provided with an electric gate or an overflow weir are sequentially arranged in an open channel; the water pump group is connected with the detection controller, and the detection controller is provided with a flow measurement unit, a water level measurement unit, a flow velocity measurement unit, a wave monitoring unit, an edge calculation unit, a communication network module, a clock synchronization unit and an execution control module; the dynamic stability of the open channel-pump coupling system is quantified by introducing a Lyapunov index and an attractor stability index (ASI), constraints such as stability, Froude number, flow velocity uniformity and the like are incorporated into an optimization framework, and the unified targets of verification precision guarantee, energy consumption optimization and hydraulic jump prevention are achieved.

Description

Technical Field

[0001] This invention relates to an open channel flow rate determination system and a multi-pump collaborative method based on attractor stability, which is a hydraulic system and control method, and a detection and regulation system and method for water conveyance and irrigation channels. Background Technology

[0002] Open channel flow metering plays a crucial role in water resource management, agricultural irrigation, and hydrological monitoring. Open channel flow meter calibration devices are key facilities for ensuring the accuracy of flow measurement, and their core requirement is to establish a stable and uniform flow field in the measurement section. The core shortcomings of existing technologies are:

[0003] There is a lack of dynamic stability assessment methods for open channel-pump coupled systems: existing controllers have not established a dynamic model that couples the open channel free water surface, pump characteristics, and gate characteristics, making it impossible to quantitatively assess system stability.

[0004] Froude number constraint not considered: Traditional methods only control flow rate or water level, without taking Fr as an optimization constraint, which may cause the system to operate near the critical flow, which is extremely unstable.

[0005] Passive response to water jumps and waves: Existing methods only take measures after detecting drastic fluctuations in water level, which is a reactive measure and may have already disrupted the verification process.

[0006] Active control that ignores velocity uniformity: Velocity uniformity mainly relies on passive rectification devices (guide grids, honeycomb grids) in the steady flow section, and is not used as an optimization target for pump control.

[0007] The decision on the number of pumps lacks a theoretical basis: when to start and stop the pumps depends on empirical thresholds and does not consider the impact of different pump configurations on the dynamic stability of the system.

[0008] Therefore, there is an urgent need for a new control method that can establish a dynamic model and stability evaluation index for the open channel-pump coupled system, actively optimize the Froude number, prevent hydraulic jumps, predictively avoid surface wave interference, incorporate flow velocity uniformity into the optimization objective, and achieve global coordinated optimization of pump number, frequency, and gate opening. Summary of the Invention

[0009] To overcome the problems of existing technologies, this invention proposes an open channel flow rate calibration system and a multi-pump collaborative method based on attractor stability. Addressing issues such as flow field instability, hydraulic jump risk, and high energy consumption during open channel flowmeter calibration, this invention establishes a dynamic model of the open channel-multi-pump coupled system. It innovatively introduces the Attractor Stability Index (ASI) to quantify system stability, incorporating open channel-specific constraints such as Froude number and flow velocity uniformity into the optimization framework, achieving global collaborative optimization of calibration accuracy, energy consumption, and stability. The system and method described in this invention quantify the dynamic stability of the open channel-pump coupled system by introducing the attractor stability index, incorporating stability constraints and open channel-specific parameters such as flow velocity uniformity and Froude number into the optimization objectives. This achieves optimization of energy consumption and system stability while ensuring calibration accuracy, making it particularly suitable for open channel flowmeter calibration devices, irrigation canal flow control, and hydrological station measurements.

[0010] The objective of this invention is achieved as follows: a flow rate determination system for open channels based on attractor stability, comprising: a water storage tank, a pump group with at least two pumps connected in parallel, a flow straightening section equipped with a guide grid or flow straightening grid, a flow stabilization section, a measurement section, and an outlet section equipped with an electric gate or overflow weir, all arranged sequentially in the open channel; the pump group is connected to a detection controller, which includes a flow measurement unit, a water level measurement unit, a flow velocity measurement unit, a wave monitoring unit, an edge computing unit, a communication network module, a clock synchronization unit, and an execution control module; the detection controller is connected to a cloud server.

[0011] The flow measurement unit is installed in the outlet pipe of each pump;

[0012] The water level measurement unit is installed at one location each upstream, middle, and downstream of the measurement section;

[0013] The flow velocity measurement unit is installed at points 5-9 on the cross section of the measurement section;

[0014] The edge computing unit includes:

[0015] A Fakubi matrix calculation module using numerical differentiation or automatic differentiation is employed.

[0016] An eigenvalue solving module using the QR decomposition algorithm or the Jacobi iterative algorithm is employed.

[0017] A multi-objective optimization solution module employing genetic algorithm, particle swarm optimization algorithm, or NSGA-II algorithm is used.

[0018] A simulation module for unsteady open channel flow based on the finite difference method or finite volume method of the Saint-Venant equations;

[0019] A stability map storage module that stores pre-computed M(n, Qd) data in the form of a database or memory-mapped file.

[0020] The cloud server provides historical data storage services for storing ASI time series, Froude number time series, flow rate distribution, and calibration results; it also provides stability map update services, optimizes M parameters based on long-term operating data, and automatically generates calibration reports.

[0021] The remote monitoring visualization interface displays the channel longitudinal profile water level line, flow velocity vector field, and attractor trajectory.

[0022] Furthermore, the edge computing unit includes:

[0023] The fast Lyapunov exponent estimation module employs the following algorithm:

[0024] Extract the most recent M state vectors {x(t1), x(t2), ..., x(tM)} from the historical buffer, where M = 80-150;

[0025] Calculate the distance between adjacent states:

[0026]

[0027] Normalization of variables: Flow normalization: Water level normalization: Flow velocity normalization: ;

[0028] Calculate the normalized distance growth rate:

[0029]

[0030] This estimated value serves as a real-time stability monitoring indicator.

[0031] ASI Fast Approximation Module:

[0032] Calculate the Jacobian matrix at the current state x0;

[0033] The Arnoldi iteration only calculates the 4-6 eigenvalues ​​with the largest modulus;

[0034] Use the formula:

[0035]

[0036] Calculation time: 25-40ms;

[0037] Module for rapid assessment of flow velocity uniformity;

[0038] The cross-sectional velocity distribution {v1, v2, ..., vm} is measured using the 5-point method or the 9-point method.

[0039] Calculate the flow velocity non-uniformity:

[0040]

[0041] Real-time determination of whether it meets the requirements .

[0042] A multi-pump coordinated method using the above-mentioned open channel flow measurement system based on attractor stability, the method comprising the following steps:

[0043] Step 1, Dynamic Modeling of the Open Channel-Pump System: Based on the measured characteristic curves of each pump and the hydraulic characteristics of the open channel, establish a set of state-space dynamic equations:

[0044]

[0045] Where: Qi is the instantaneous flow rate of the i-th pump; Where: is the pump hydraulic time constant, in seconds; Qtarget,i is the target flow rate determined by the frequency fi and the channel water level h in the calibration section; kc is the pump coupling coefficient; Ac is the channel cross-sectional area; Qout is the outflow flow rate, determined by the downstream gate opening and water level; v is the average flow velocity in the measurement section; L is the length of the steady flow section; g is the acceleration due to gravity; h0 is the channel bottom elevation; n is the Manning roughness coefficient; R is the hydraulic radius; N is the number of currently operating pumps.

[0046] Step 2, real-time evaluation of attractor stability:

[0047] Calculate the Jacobian matrix at the current running state x0 = [Q1, Q2, ..., QN, h, v]T:

[0048]

[0049] Each element is calculated using numerical differentiation:

[0050]

[0051] The step size δ is taken as the value of the corresponding variable. times;

[0052] Calculate all eigenvalues ​​{λ1, λ2, ..., λn} of the Jacobian matrix;

[0053] Define the attractor stability index ASI:

[0054]

[0055] Where β is the weighting coefficient, with a value ranging from 0.05 to 0.15, and an optimal value of 0.10; ε is the regularization constant, with a value of 0.01;

[0056] Step 3, Dynamic optimization of pump combination number:

[0057] Query the pre-built stability map M(n, Qd), which stores the ASI values ​​for different combinations of pump number n and calibration flow rate Qd;

[0058] Based on the current calibration flow rate requirement Qd, select the number of pumps that maximizes ASI.

[0059]

[0060] Pump start / stop switching is performed when any of the following conditions are met:

[0061] The current ASI is < 0.12;

[0062] The expected energy savings after the switch will exceed 8%;

[0063] |n* - ncurrent| ≥ 1 and the switching has been running stably for more than 300 seconds;

[0064] Froude number Exceeding the safe range [0.3, 0.7];

[0065] Step 4, Frequency Allocation Multi-Objective Optimization:

[0066] After determining the number of pumps n*, solve the optimization problem:

[0067]

[0068] Constraints include: Flow balance constraints: Frequency range constraint: fmin ≤ fi ≤ fmax; Water level range constraint: hmin ≤ h ≤ hmax; Flow velocity uniformity constraint: Stability constraints: Its resonance avoidance constraint: For all identified resonant frequencies fres,k; Froude number constraint: 0.3 < Fr < 0.7; frequency discretization constraint: Step length The preferred frequency is 1Hz; the weighting coefficients w1, w2, and w3 are set according to the operating mode: Energy-saving mode: w1=0.6, w2=0.3, w3=0.1; Precision mode: w1=0.2, w2=0.3, w3=0.5.

[0069] Step 5, Pre-verification of hydraulic transient risks:

[0070] Input the current state and the optimized target frequency f* into the open channel unsteady flow calculation model;

[0071] Solve using the Saint-Venant equations:

[0072]

[0073] Calculate the water level change amplitude Δhmax and the flow velocity change rate dv / dt during the flow regulation process:

[0074] Assessing the risk of a hydraulic jump:

[0075] If dFr / dt > 0.05 / s or the predicted Fr > 0.8, a hydraulic jump risk is identified, and the following measures should be taken:

[0076] Extend the adjustment time to ;

[0077] Prioritize adjusting the downstream gate opening, and gradually adjust it in conjunction with the pump frequency;

[0078] Phased adjustment: First adjust one pump, and adjust the other pumps after it stabilizes;

[0079] Determining wave propagation:

[0080] Calculate the surface wave propagation speed If the flow rate adjustment rate dv / dt > c / L, it may cause wave interference, and the adjustment should be delayed.

[0081] Step 6, Control command execution and feedback monitoring:

[0082] Send frequency command fi* and corresponding ramp time Δt to each inverter:

[0083] Synchronous adjustment of downstream gate opening To maintain water level stability;

[0084] Real-time acquisition of data from flow rate, water level, and flow velocity sensors, with a sampling period of no more than 1 second;

[0085] Add the newly collected data to the historical buffer, with the buffer length set to 50-200 sampling points;

[0086] Calculate the short-term Lyapunov exponent estimate within the sliding window:

[0087]

[0088] Where: d(t) is the distance between adjacent trajectory points in phase space, and Twindow is the window duration, which is 5-10 times. ;

[0089] When continuously detected Emergency stabilization control is triggered if the duration exceeds 5 seconds, or if the following anomalies are detected:

[0090] Water level fluctuations exceeding ±5cm and frequency >0.5Hz;

[0091] The Froude number has exceeded 0.75;

[0092] Flow velocity non-uniformity in the measurement section >5%;

[0093] Emergency control measures:

[0094] Freeze all frequency adjustment commands immediately:

[0095] Query the historical buffer to find the most recent one that satisfies the condition. The state fsafe;

[0096] Asymptotically recover to fsafe with maximum ramp time;

[0097] Synchronously adjust the downstream gate to the corresponding opening degree;

[0098] Send an alarm signal to the monitoring system;

[0099] Return to step 2 and continue executing the loop at 1-second intervals until one of the following stopping conditions is met, at which point the loop stops:

[0100] The current testing site has maintained stable compliance for a longer period than the testing sampling time.

[0101] The operator issues a stop command through the monitoring interface or the local control cabinet;

[0102] Single-point stability timeout;

[0103] An emergency fault has been detected.

[0104] Furthermore, the stability map M(n, Qd) is generated through the following offline calculation steps:

[0105] A1, Parameter space discretization:

[0106] Define the computational grid:

[0107] Pump count range: ;

[0108] Verification flow range: The sampling interval is segmented according to the measurement range:

[0109] Low flow rate range (0-50 L / s): 5 L / s interval;

[0110] Medium flow rate range (50-500 L / s): 20 L / s interval;

[0111] High flow rate range (500-2000 L / s): 50 L / s interval;

[0112] For each number of pumps n, generate a frequency combination set Fn, requiring f1 ≥ f2 ≥ ... ≥ fn, with a frequency step size of 2Hz;

[0113] A2. Steady-state hydraulic calculation of open channels:

[0114] For each parameter combination (n, Qd, f):

[0115] Calculate steady-state water level using the formula for uniform flow in an open channel:

[0116] ;

[0117] Where: S is the slope of the canal bottom;

[0118] Obtain the operating points (Qi, Hi) of each pump based on the pump characteristic curve;

[0119] The velocity distribution in the measurement section is calculated using either the logarithmic law or the power law.

[0120] Verify whether the Froude number Fr is within a safe range;

[0121] A3. Linearization of dynamics:

[0122] In steady-state solution Construct the Jacobian matrix at:

[0123] ;

[0124] Using the central difference scheme:

[0125]

[0126] Where: ej is the j-th standard basis vector, and the step size is... Take the scale of the variable times;

[0127] For open channel systems, the key partial derivatives include:

[0128] The water level-velocity coupling term reflects the effect of gravity. The downstream outflow characteristics depend on the gate type; : Pump operating point drift with water level changes

[0129] A4. Stability Quantification Calculation:

[0130] Calculate the eigenvalue spectrum {λ1, λ2, ..., λn} of the Jacobian matrix;

[0131] Calculate the ASI value;

[0132] Constraint violation conditions: If the water level exceeds [hmin, hmax], it is marked as infeasible; if the Froude number Fr < 0.2 or Fr > 0.8, it is marked as a risky condition; if the flow velocity non-uniformity in the measurement section is > 3%, it is marked as low accuracy; if ASI < 0, it is marked as unstable; if the efficiency of any pump is < 40%, it is marked as inefficient.

[0133] A5. Data Storage and Indexing:

[0134] Calculation results Save to an SQLite database or an HDF5 file;

[0135] Create a multidimensional index: Primary key: (n, Qd); Indexed fields: .

[0136] 5. The method according to claim 4, wherein the resonant frequency fres is obtained by the following online identification method:

[0137] B1. Small-amplitude incentive injection:

[0138] When the system is running stably (selecting a medium flow rate condition, Fr≈0.5), inject a sinusoidal frequency disturbance into the selected pump:

[0139] fi(t) = fi,0 + A·sin(2πftest·t)

[0140] Wherein: the reference frequency fi,0 is the current operating frequency, the amplitude A is 0.3-0.8Hz, the test frequency ftest is scanned from 0.05Hz to 5Hz, and the duration of each test frequency is 30-60 seconds;

[0141] B2. System Response Measurement:

[0142] Simultaneously acquire the following at a sampling rate of not less than 20Hz: water level signal hw(t) in the middle of the measurement section; average flow velocity signal v(t) in the measurement section; optional: water surface fluctuation signal. ;

[0143] For each test frequency, extract the steady-state oscillation data and calculate the frequency response function:

[0144] ;

[0145] B3. Resonance Peak Identification:

[0146] Detect frequency points fr that meet one of the following conditions in the frequency sweep results:

[0147] Water level resonance criterion: The amplitude is a local maximum, and Half-power bandwidth < 0.2Hz (low Q value due to open channel resonance); phase jump > 45° near fr;

[0148] Velocity resonance criterion: ; accompanied by increased water surface ripples;

[0149] The identified fr frequency is recorded as the system's natural frequency. Typical resonance sources include: channel standing wave resonance. ; Inlet pool sloshing: fr = 0.1-0.5 Hz; Pump blade passing frequency harmonics;

[0150] B4. Frequency zone setting is disabled:

[0151] For each resonant frequency fres,k, establish a forbidden interval:

[0152] ;

[0153] Among them, safety boundary Select 1.5-2.5Hz;

[0154] The optimization constraints in step S4 are as follows: the pump frequency fi must not fall within Fforbidden; the pump frequency difference... The frequency must not fall into Fforbidden; the 2nd and 3rd harmonics of the pump frequency must not fall into Fforbidden.

[0155] Furthermore, the target flow rate Qtarget,i(fi, h) is calculated based on the following pump characteristic model and open channel hydraulic coupling:

[0156] Pump characteristic curve:

[0157]

[0158] The optimal efficiency point flow rate (QBEP) varies linearly with frequency.

[0159] QBEP(fi) = QBEP,rated·(fi / frated);

[0160] Open channel hydraulic coupling:

[0161] The pump head needs to overcome: channel water level h; pipeline friction loss. ;

[0162] Local loss ;

[0163] Overall head balance equation:

[0164] Hi = h + hf + hm + hvp

[0165] Where: hvp is the difference between the pump installation elevation and the canal bottom elevation;

[0166] The target flow Qtarget,i is the solution to the following system of equations:

[0167] .

[0168] Furthermore, the following adaptive parameter update steps are also included:

[0169] Periodic parameter calibration:

[0170] Collect operational log data, including frequency commands, measured flow rate, measured water level, and measured flow velocity distribution;

[0171] Update dynamic model parameters using system identification algorithms :

[0172] Use either an extended Kalman filter (EKF) or an unscented Kalman filter (UKF);

[0173] Fit the measured water level response curve to the model prediction residual;

[0174] Update Manning roughness coefficient n:

[0175] Extract (Qsteady, hsteady, vsteady) from steady-state operating data;

[0176] Using Manning's formula to calculate inversely: ;

[0177] Detect n-value drift;

[0178] Update the pump characteristic curve parameters H0,rated, KH,rated, ηmax:

[0179] Extract steady-state operating points (Qi, Hi, fi) from the operating data;

[0180] The characteristic curve equation was fitted using multivariate nonlinear regression.

[0181] Inspect impeller wear;

[0182] Update stability map M:

[0183] Recalculate the ASI values ​​for the critical region using the new parameters;

[0184] An incremental update strategy is adopted, and only areas where changes exceed 15% are recalculated.

[0185] The original map is retained as a backup, and the replacement is performed gradually.

[0186] Online learning optimization:

[0187] Record each verification task (Qd, n, f) and the verification result:

[0188] Establish a decision-outcome database;

[0189] Gradual optimization using reinforcement learning:

[0190] State space: ;

[0191] Action space: (n, f1, f2, ..., fn);

[0192] Reward function: R = -uncertainty - 0.1·energy consumption + 10·(ASI>0.15);

[0193] Online correction of channel hydraulic characteristics:

[0194] Detecting siltation at the bottom of the canal:

[0195] If the water level rises systematically under the same flow rate, it is inferred that the cross-sectional area of ​​the water passage has decreased.

[0196] Corrected effective cross-sectional area ;

[0197] Recalculate the Froude number threshold;

[0198] Detect changes in gate characteristics:

[0199] Identify the gate flow coefficient Cd and fit it to the measured Qout-h relationship; detect gate corrosion or deformation.

[0200] Furthermore, the calculation of unsteady flow in the open channel employs the following method:

[0201] Simplified calculation based on the method of characteristics:

[0202] For rectangular channels, a single-wave approximation is used:

[0203] ;

[0204] Furthermore, the calculation of unsteady flow in the open channel employs the following method:

[0205] Detailed simulation based on Saint-Venant's equations:

[0206] Discretize the Saint-Venant equations using the Preissmann four-point implicit difference scheme:

[0207] :

[0208]

[0209] Boundary conditions: Upstream: Downstream: The hQ relationship is determined by the gate equation;

[0210] Calculation output: water level time history h(x,t) for each cross-section; flow velocity time history v(x,t) for each cross-section; Froude number time history Fr(x,t); flow velocity nonuniformity time history σv(t) for the measurement section;

[0211] Risk assessment: If max[Fr(x,t)] > 0.85, there is a risk of hydraulic jump; if The water level changes too quickly; if The flow field stability is insufficient.

[0212] Furthermore, the optimization problem is solved using the following algorithm:

[0213] Constrained particle swarm optimization;

[0214] Particle definition: Each particle represents a frequency combination f = [f1, f2, ..., fn];

[0215] Speed ​​updates:

[0216]

[0217] Location update:

[0218]

[0219] Constraint handling: The penalty function method is used, and the augmented objective function is defined as follows:

[0220]

[0221] Where: M is the large penalty factor, and its value is... ;

[0222] Parameter settings;

[0223] Number of particles: 30-50; Number of iterations: 40-60 generations; Inertia weight. ;

[0224] Learning factors: c1 = c2 = 2.0.

[0225] Furthermore, the optimization problem is solved using the following algorithm:

[0226] Multi-objective genetic algorithm method:

[0227] When it is necessary to optimize energy consumption, accuracy, and stability simultaneously, NSGA-II should be used.

[0228] Encoding: Integer encoding, gene = frequency rank number;

[0229] Fitness: Multi-objective vector ;

[0230] Pareto sort: Non-dominated sort + crowding distance;

[0231] Genetic operators: Selection: Tournament selection; Crossover: Simulated binary crossover, probability 0.9; Mutation: Polynomial mutation, probability 0.1; Output: Multiple solutions on the Pareto front.

[0232] Furthermore, the optimization problem is solved using the following algorithm:

[0233] Sequential quadratic programming method:

[0234] Suitable for small-scale problems with ≤3 pumps, with fast convergence speed;

[0235] Solve the quadratic programming subproblem in each iteration:

[0236] ;

[0237] Update the Hessian matrix approximation using the BFGS formula;

[0238] Optimal algorithms: NSGA-II for precision mode, PSO for energy-saving mode, and SQP for real-time response.

[0239] The advantages and beneficial effects of this invention are as follows: This invention quantifies the dynamic stability of the open channel-pump coupling system by introducing the Lyapunov index and the attractor stability index (ASI), and incorporates constraints such as stability, Froude number, and flow velocity uniformity into the optimization framework, thereby achieving the unified goals of ensuring calibration accuracy, optimizing energy consumption, and preventing hydraulic jumps.

[0240] This invention is the first to apply dynamic system theory to an open channel flow calibration system and introduces the ASI index to quantify the stability of the open channel-pump coupling system, filling a theoretical gap in this field. Traditional methods rely on empirical parameters (such as "water level fluctuation < 5 cm"), while ASI provides a numerically defined index with clear physical meaning. This invention proactively prevents hydraulic jumps, ensuring a high calibration success rate. It predicts Fr changes using the Saint-Venant equation and adjusts the control strategy before Fr approaches 1.0, preventing hydraulic jumps from disrupting the flow field in the measurement section. Traditional methods can only detect hydraulic jumps retrospectively (through a sudden rise in water level), by which time the calibration has already failed.

[0241] This invention incorporates flow velocity uniformity into the optimization objective, relying not only on passive rectification in the steady flow section but also on actively improving the flow velocity distribution by optimizing pump frequency combinations. Theoretical analysis shows that different frequency combinations produce different pulsation superposition effects, and optimization can reduce flow velocity non-uniformity from 3-5% to <2%.

[0242] This invention globally and collaboratively optimizes the number of pumps, frequency, and gate opening, combining the three for optimization. Compared with traditional hierarchical control (first determining the number of pumps, then adjusting the frequency, and finally adjusting the gate), it can reduce energy consumption by 8-15% while improving stability.

[0243] This invention adapts to the unique physical constraints of open channels and considers constraints that do not exist in pipe systems: Froude number range (avoiding critical flow and rapid flow); surface wave propagation (regulating velocity limits); channel bottom friction (Manning formula); free surface fluctuations (more sensitive than pipes).

[0244] This invention has strong parameter adaptability and can identify Manning roughness coefficient n (to detect siltation or algae adhesion), gate flow coefficient Cd (to detect corrosion or deformation), and pump characteristic parameters (to detect impeller wear) online.

[0245] This invention is theoretically reliable and engineeringally feasible. Based on the classical Saint-Venant equations and dynamical system theory, it has a solid theoretical foundation. The required hardware (ADCP, level gauge, frequency converter) are all mature products, and implementation costs are controllable. Theoretical simulations (based on HEC-RAS software) show that compared to traditional PID control, flow velocity uniformity is improved by 30-50%, the probability of hydraulic jump is reduced from 15% to <1%, the stabilization time is shortened from 5-8 minutes to 2-3 minutes, and energy consumption is reduced by 10-18%. Attached Figure Description

[0246] The present invention will be further described below with reference to the accompanying drawings and embodiments.

[0247] Figure 1 This is a schematic diagram of the open channel flow rate verification system for attractor stability as described in Embodiment 1 of the present invention.

[0248] Figure 2 This is a flowchart of the method described in Embodiment 2 of the present invention. Detailed Implementation

[0249] Example 1:

[0250] This embodiment is an open channel flow measurement system based on attractor stability, such as... Figure 1 As shown.

[0251] Open channel flow metering plays a vital role in water resource management, agricultural irrigation, and hydrological monitoring. Open channel flow meter calibration devices are crucial facilities for ensuring the accuracy of flow measurement; their core requirement is to establish a stable and uniform flow field in the measurement section.

[0252] Typical structure of open channel flow measurement device:

[0253] Water storage tank → Multiple pumps in parallel supply → Rectifying section → Flow stabilizing section → Measuring section (installation of the flow meter under test) → Outlet section (gate / overflow weir).

[0254] The requirements for the calibration device are: the uniformity of the flow velocity distribution in the measurement section is less than 2%, the flow fluctuation is less than 1%, and the uncertainty is less than or equal to 0.5%. Based on these flow measurement requirements, this embodiment adds various sensors and computing facilities, including a data acquisition module (flowmeter, water level gauge, and flow profiler (ADCP)); an edge computing unit (integrated dynamic model, optimization algorithm, and Saint-Venant equation solver); and an execution control module (frequency converter and gate actuator). A dynamic model of the open channel-multi-pump coupled system is established, considering the coupling of free surface, Froude number, and Manning friction. This results in an open channel flow calibration device based on the attractor stability of the dynamic system, thus forming a multi-pump collaborative control system.

[0255] The open channel flow calibration system described in this embodiment includes: a data acquisition module (flow meter, water level gauge, and flow velocity profiler (ADCP)); an edge computing unit (integrated dynamic model, optimization algorithm, and Saint-Venant equation solver); and an execution control module (frequency converter and gate electric actuator). The open channel is sequentially configured with a reservoir, a pump set with at least two parallel pumps, a channel straightening section equipped with a guide grid or flow straightening grid, a flow stabilization section, a measurement section, and an outlet section equipped with an electric gate or overflow weir. The pump set is connected to a detection controller, which includes a flow measurement unit, a water level measurement unit, a flow velocity measurement unit, a wave monitoring unit, an edge computing unit, a communication network module, a clock synchronization unit, and an execution control module. The detection controller is connected to a cloud server via a wired or wireless network.

[0256] The flow measurement unit uses an electromagnetic flow meter or an ultrasonic flow meter with a measurement accuracy better than 0.3%FS and is installed in the outlet pipe of each pump.

[0257] The water level measurement unit adopts an ultrasonic water level gauge or a pressure water level gauge with a measurement accuracy better than ±1mm, and is installed at one location each upstream, middle and downstream of the measurement section.

[0258] The flow velocity measurement unit employs an ultrasonic Doppler velocity profiler or a point velocity meter, with a measurement accuracy better than ±0.5%, and is deployed at 5-9 points on the cross-section of the measurement section. The wave monitoring unit uses a high-frequency ultrasonic sensor with a sampling rate ≥10Hz. The edge computing unit includes an embedded processor with a main frequency of not less than 1.5GHz, memory of not less than 1GB, and running a real-time operating system. It also includes a Fakubi matrix calculation module using numerical differentiation or automatic differentiation; an eigenvalue solving module using QR decomposition or Jacobi iterative algorithms; a multi-objective optimization solving module using genetic algorithms, particle swarm optimization algorithms, or NSGA-II algorithms; a finite difference method or finite volume method for open channel unsteady flow simulation based on the Saint-Venant equations; and a stability map storage module that stores pre-calculated M(n, Qd) data in the form of a database or memory-mapped file.

[0259] The communication network module includes: a fieldbus using Modbus RTU / TCP, Profinet, or EtherCAT protocols; and an optional wireless communication module using 4G / 5G or LoRa technology. A clock synchronization unit using NTP or PTP protocols ensures multi-point measurement timestamp accuracy better than 5ms.

[0260] The execution control module includes: a frequency converter with a frequency resolution better than 0.1Hz and a frequency response time of less than 200ms; a pump start / stop control unit using a soft starter or electromagnetic contactor; a downstream gate electric actuator with an opening resolution better than 0.5° and a response time of <5 seconds; and a local control cabinet supporting manual / automatic / calibration modes.

[0261] The cloud server provides historical data storage services for storing ASI time series, Froude number time series, flow velocity distribution, and verification results; it provides stability map update services, optimizes M parameters based on long-term operating data, and automatically generates verification reports in accordance with the JJG 164-2000 standard format; it provides a remote monitoring visualization interface to display channel longitudinal profile water level lines, flow velocity vector fields, and attractor trajectories.

[0262] The edge computing unit further includes:

[0263] The fast Lyapunov exponent estimation module employs the following algorithm:

[0264] Extract the most recent M state vectors {x(t1), x(t2), ..., x(tM)} from the historical buffer, where M = 80-150 (open channel systems are dynamic and require a longer window).

[0265] Calculate the distance between adjacent states:

[0266] .

[0267] Normalization of variables: Flow normalization: Water level normalization: Flow velocity normalization: .

[0268] Calculate the normalized distance growth rate:

[0269]

[0270] This estimate serves as a real-time stability monitoring indicator.

[0271] ASI Fast Approximation Module:

[0272] Calculate the Jacobian matrix at the current state x0:

[0273] The Arnoldi iteration only calculates the 4-6 eigenvalues ​​with the largest modulus.

[0274] Using an approximate formula:

[0275] .

[0276] Computation time: 25-40ms (slower than pipeline systems because the state dimension includes h and v).

[0277] Rapid assessment module for flow rate uniformity:

[0278] The cross-sectional velocity distribution {v1, v2, ..., vm} is measured using the 5-point method or the 9-point method.

[0279] Calculate the flow velocity non-uniformity:

[0280]

[0281] Real-time determination of compliance with ISO 4359 standard ( ).

[0282] Example 2:

[0283] This embodiment describes a multi-pump collaborative method for an open channel flow rate verification system using the attractor stability described in the above embodiments.

[0284] The method employs the following technical solution:

[0285] Firstly, a multi-pump collaborative control method for open channel flow calibration device based on the stability of attractors in dynamic systems is provided. The core idea is to model the open channel-pump coupled system as a dynamic system in state space, evaluate and ensure stability by analyzing the attractor characteristics of the system, and simultaneously satisfy the Froude number constraint and flow velocity uniformity requirements of open channel flow.

[0286] Specifically, the method comprises six main steps:

[0287] (1) Dynamic modeling of open channel-pump system: Establish a state-space equation set including pump flow rate, channel water level, and flow velocity in the measurement section. The model considers: hydraulic coupling between pumps (shared intake pool), coupling between pump and channel (pump head needs to overcome water level), unsteady flow characteristics of the channel (interaction between water level, flow rate, and flow velocity), and feedback effect of downstream gate.

[0288] (2) Attractor stability assessment: Linearize the dynamic equations at the current operating point, calculate the eigenvalue spectrum of the Jacobian matrix, and define the attractor stability index (ASI). Compared to pipeline systems, the ASI calculation for open channel systems requires the inclusion of a coupling term for water level h and flow velocity v.

[0289] (3) Dynamic optimization of pump number: Based on the calibration flow requirements, query the pre-established stability map and select the pump number configuration with the maximum ASI and satisfying the Fr constraint.

[0290] (4) Frequency allocation multi-objective optimization: Establish a multi-objective optimization problem that includes energy consumption, flow deviation and flow velocity uniformity, with constraints including ASI threshold, Fr range, resonance avoidance, etc.

[0291] (5) Pre-verification of unsteady flow in open channels: The Saint-Venant equations are used to predict water level changes and Froude number changes during flow regulation. If a hydraulic jump risk is predicted (Fr close to 1.0), the regulation time is automatically extended or the gate opening is adjusted.

[0292] (6) Execution and monitoring: After issuing control commands, monitor ASI, Fr, and flow velocity distribution in real time. If deviation from the stable attractor or Fr abnormality is detected, trigger emergency stabilization control.

[0293] Secondly, an offline method for generating stability maps is provided, taking into account the characteristics of open channel systems: the steady state is calculated using the uniform flow formula for open channels; the Fr at each operating point is verified to be within the safe range; and the uniformity of flow velocity is calculated (using logarithmic law or power law).

[0294] Thirdly, an online method for identifying resonant frequencies is provided, taking into account the characteristics of open channels: monitoring both water level and flow velocity responses; identifying the resonant frequencies of standing waves in the channel; and ensuring that the restricted area is wider than that of a pipeline system (due to the low damping of open channels).

[0295] This embodiment primarily employs a five-step control strategy: Offline mapping: Pre-calculating ASI values ​​under different combinations of pump numbers, flow rates, and frequencies to generate a stability map. Pump number decision: Based on the verification flow rate requirements, querying the map to select the pump combination with the highest ASI that satisfies the Fr constraint. Frequency optimization: Establishing a multi-objective optimization model (energy consumption + flow deviation + velocity uniformity), with constraints including ASI threshold, Froude number range, Manning friction balance, and resonance avoidance. Hydraulic pre-verification: Using the Saint-Venant equation to simulate the regulation process and predict the Fr trajectory; if a hydraulic jump risk is detected, extending the regulation time or coordinating gate operation. Real-time monitoring: Calculating the short-term Lyapunov exponent and real-time ASI to detect early signs of system deviation from the stable attractor; triggering emergency backoff when an anomaly occurs.

[0296] The specific steps of the method are as follows:

[0297] Step 1, Dynamic modeling of the open channel-pump system:

[0298] Based on the measured characteristic curves of each pump and the hydraulic characteristics of the open channel, a set of state-space dynamic equations is established:

[0299]

[0300] Where: Qi is the instantaneous flow rate of the i-th pump, in m³ / s; 1. Pump hydraulic time constant, in seconds, determined by step response test, with a value ranging from 0.5 to 3 seconds; 2. Target,i, the target flow rate determined by frequency fi and channel water level h; 3. kc, the pump coupling coefficient, ranging from 0.01 to 0.1; 4. h, the water level in the calibration section, in meters; 5. Ac, the channel cross-sectional area, in m²; 6. Qout, the outflow rate, determined by the downstream gate opening and water level; 7. v, the average flow velocity in the measurement section, in m / s; 8. L, the length of the steady flow section, in meters; 9.81 m / s², the gravitational acceleration; 10. h0, the channel bottom elevation, in meters; 11. n, the Manning roughness coefficient, ranging from 0.012 to 0.025; 2. R, the hydraulic radius, in meters; 3. N, the number of currently operating pumps.

[0301] Step 2, real-time evaluation of attractor stability:

[0302] Calculate the Jacobian matrix at the current running state x0 = [Q1, Q2, ..., QN, h, v]T:

[0303]

[0304] Each element is calculated using numerical differentiation:

[0305]

[0306] The step size δ is taken as the value of the corresponding variable. times.

[0307] Calculate all eigenvalues ​​of the Jacobian matrix. ;

[0308] Define the attractor stability index ASI:

[0309]

[0310] Wherein: β is the weighting coefficient, with a value range of 0.05-0.15 and an optimal value of 0.10; This is the regularization constant, with a value of 0.01. Let i be the eigenvalues, i = 1 to n.

[0311] Step 3, Dynamic optimization of pump combination number:

[0312] Query the pre-built stability map M(n, Qd), which stores the ASI values ​​for different combinations of pump number n and calibration flow rate Qd.

[0313] Based on the current calibration flow rate requirement Qd, select the number of pumps that maximizes ASI:

[0314]

[0315] Pump start / stop switching is performed when any of the following conditions are met: current ASI < 0.12, expected energy savings exceeding 8% after switching, |n* - ncurrent| ≥ 1 and the switching has been running stably for more than 300 seconds, and Froude number. Exceeding the safe range [0.3, 0.7].

[0316] Step 4, Frequency Allocation Multi-Objective Optimization:

[0317] After determining the number of pumps n*, solve the optimization problem:

[0318]

[0319] Constraints include: Flow balance constraints: Frequency range constraint: fmin ≤ fi ≤ fmax; Water level range constraint: hmin ≤ h ≤ hmax; Flow velocity uniformity constraint: (Measurement section velocity fluctuation coefficient), stability constraints: ,in Values ​​range from 0.10 to 0.18, with a preferred value of 0.15; resonance avoidance constraint: For all identified resonant frequencies fres,k, the Froude number constraint is 0.3 < Fr < 0.7 (to avoid critical flow and rapid flow), and the frequency discretization constraint is: Step length 1Hz is preferred.

[0320] The weighting coefficients w1, w2, and w3 are set according to the operating mode: Energy saving mode: w1=0.6, w2=0.3, w3=0.1; Precision mode: w1=0.2, w2=0.3, w3=0.5.

[0321] Step 5, Pre-verification of hydraulic transient risks:

[0322] Input the current state and the optimized target frequency f* into the open channel unsteady flow calculation model;

[0323] Solve using the Saint-Venant equations:

[0324]

[0325] Calculate the water level change amplitude during flow regulation. And the rate of change of flow velocity dv / dt.

[0326] Assessing the risk of a hydraulic jump:

[0327] If dFr / dt > 0.05 / s or the predicted Fr > 0.8, a hydraulic jump risk is identified, and the following measures are implemented: extend the slope adjustment time to Prioritize adjusting the downstream gate opening, and adjust the pump frequency slowly in conjunction with the adjustment; adjust in stages: adjust one pump first, and adjust the other pumps after stabilization.

[0328] Determining wave propagation:

[0329] Calculate the surface wave propagation speed If the flow rate adjustment rate dv / dt > c / L, it may cause wave interference, and the adjustment needs to be delayed.

[0330] Step 6, Control command execution and feedback monitoring:

[0331] Send frequency command fi* and corresponding ramp time to each inverter Synchronously adjust the opening of the downstream gate. To maintain stable water levels, real-time data from flow rate, water level, and flow velocity sensors are collected, with a sampling period of no more than 1 second. Newly collected data is added to the historical buffer, with the buffer length set to 50-200 sampling points.

[0332] Calculate the short-term Lyapunov exponent estimate within the sliding window:

[0333]

[0334] Where: d(t) is the distance between adjacent trajectory points in phase space, and Twindow is the window duration, which is 5-10 times. .

[0335] When continuously detected Emergency stabilization control is triggered if the duration exceeds 5 seconds, or if the following anomalies are detected:

[0336] Water level fluctuations exceeding ±5cm and frequency >0.5Hz (signs of waves), Froude number exceeding 0.75 (approaching critical flow), and flow velocity non-uniformity in the measurement section >5% (intensified turbulence).

[0337] Emergency control measures:

[0338] Immediately freeze all frequency adjustment commands, query the historical buffer, and find the most recent command that satisfies ASI > 0.18 and The system gradually restores the gate to its fsafe state over a maximum ramp time (preferably 15 seconds), synchronously adjusts the downstream gate to the corresponding opening, and sends an alarm signal to the monitoring system.

[0339] Return to step 2 and continue executing the loop at 1-second intervals until one of the following stopping conditions is met, at which point the loop stops:

[0340] The current testing sites consistently meet the standards (ASI>0.15). The duration of flow deviation <1% exceeds the calibration sampling time;

[0341] The operator issues a stop command through the monitoring interface or the local control cabinet;

[0342] Single-point stability timeout (failed to meet standard after more than 15 minutes);

[0343] An emergency fault was detected (pump failure, water level exceeding limit, communication interruption).

[0344] The conditions for stopping include three types: automatic stop conditions, manual stop conditions, and emergency stop conditions.

[0345] Automatic stop conditions (any one of them must be met):

[0346] 1. Verification Completion Conditions: Data collection at the current verification point is completed, and all of the following conditions are met for a duration exceeding the set verification sampling time (e.g., 60-120 seconds):

[0347] .

[0348] 2. Full-range calibration completed: All calibration flow points have been measured and recorded.

[0349] 3. Operator-initiated stop: Issue a stop command through the monitoring interface or local control cabinet (used for: equipment maintenance, troubleshooting, and non-routine operation requirements).

[0350] 4. Timeout Stop: If the stabilization time of a single calibration point exceeds the preset upper limit (e.g., 15-20 minutes) and still fails to meet the standard, the system will pause and issue an alarm.

[0351] 5. Critical fault shutdown:

[0352] Any pump failure (overcurrent, overheating, excessive vibration):

[0353] The water level exceeds the safety limit (h > 0.95·hmax or h < 0.3·hdesign);

[0354] Communication interruption lasting more than 10 seconds;

[0355] Sensor data is abnormal (outside the physically reasonable range).

[0356] The stability map M(n, Qd) above is generated through the following offline calculation steps:

[0357] A1. Parameter space discretization:

[0358] Define the computational grid:

[0359] Pump count range: n {1, 2,...,Ntotal};Verification flow range: Qd [Qmin, Qmax], sampling intervals are segmented according to the measurement range:

[0360] Small flow rate range (0-50 L / s): interval 5 L / s; Medium flow rate range (50-500 L / s): interval 20 L / s; Large flow rate range (500-2000 L / s): interval 50 L / s; For each number of pumps n, generate a frequency combination set Fn, requiring f1 ≥ f2 ≥ ... ≥ fn, with a frequency step size of 2Hz.

[0361] A2. Steady-state hydraulic calculation of open channels:

[0362] For each parameter combination (n, Qd, f):

[0363] Calculate steady-state water level using the formula for uniform flow in an open channel:

[0364]

[0365] Where: S is the channel bottom slope: obtain the operating points (Qi, Hi) of each pump, calculate the flow velocity distribution of the measurement section based on the pump characteristic curve, and verify whether the Froude number Fr is within the safe range using the logarithmic law or power law.

[0366] A3. Linearization of dynamics:

[0367] Construct the Jacobian matrix at the steady-state solution xss = [Q1,ss, Q2,ss, ..., hss, vss]T:

[0368] ;

[0369] Using the central difference scheme:

[0370]

[0371] Where: ej is the j-th standard basis vector, and the step size is... Take the scale of the variable times.

[0372] For open channel systems, the key partial derivatives include:

[0373] The water level-velocity coupling term reflects the effect of gravity.

[0374] The downstream outflow characteristics depend on the gate type;

[0375] The drift of the pump operating point as the water level changes.

[0376] A4. Stability Quantification Calculation:

[0377] Calculate the eigenvalue spectrum of the Jacobian matrix Calculate the ASI value.

[0378] Constraint violation conditions: If the water level exceeds [hmin, hmax], it is marked as infeasible; if the Froude number Fr < 0.2 or Fr > 0.8, it is marked as a risky condition; if the flow velocity non-uniformity in the measurement section is > 3%, it is marked as low accuracy; if ASI < 0, it is marked as unstable; if the efficiency of any pump is < 40%, it is marked as inefficient.

[0379] A5. Data Storage and Indexing:

[0380] The calculation results (n, Qd, f, ASI, Fr, Store the total (Ptotal, flags) in an SQLite database or an HDF5 file.

[0381] Create a multidimensional index: Primary key: (n, Qd), Indexed fields: ASI, Fr, .

[0382] Optional: Use radial basis function interpolation or Kriging interpolation to improve query resolution, especially for densifying the area near commonly used flow points during verification.

[0383] The resonant frequency fres is obtained through the following online identification method:

[0384] B1. Small-amplitude incentive injection:

[0385] When the system is running stably (selecting a medium flow rate condition, Fr≈0.5), inject a sinusoidal frequency disturbance into the selected pump:

[0386]

[0387] Wherein: the reference frequency fi,0 is the current operating frequency, the amplitude A is 0.3-0.8Hz (open channel systems are more sensitive to disturbances, and the amplitude is smaller than that of pipeline systems), and the test frequency ftest is gradually scanned from 0.05Hz to 5Hz, with each test frequency lasting 30-60 seconds.

[0388] B2. System Response Measurement:

[0389] Simultaneously acquire the following data at a sampling rate of no less than 20Hz: water level signal hw(t) in the middle of the measurement section, average flow velocity signal v(t) in the measurement section, and optionally: water surface fluctuation signal. .

[0390] For each test frequency, extract the steady-state oscillation data (remove the first 10 seconds of transient data) and calculate the frequency response function:

[0391] Hh(ftest) = |FFT[hw(t)]|f=ftest / |FFT[fi(t)]|f=ftest (Water level response)

[0392] Hv(ftest) = |FFT[v(t)]|f=ftest / |FFT[fi(t)]|f=ftest (Flow rate response)

[0393] Record the amplitudes |Hh|, |Hv| and the phases ∠Hh, ∠Hv.

[0394] B3. Resonance Peak Identification:

[0395] Detect frequency points fr that meet one of the following conditions in the frequency sweep results:

[0396] Water level resonance criteria: the amplitude is a local maximum, and |Hh(fr)| > 3×mean(|Hh|), the half-power bandwidth is <0.2Hz (the Q value of open channel resonance is low), and the phase jumps >45° near fr.

[0397] Flow velocity resonance criterion: |Hv(fr)| > 2.5×mean(|Hv|), accompanied by increased water surface ripples (if a wave sensor is present).

[0398] The identified fr frequency is recorded as the system's natural frequency. Typical resonance sources include: channel standing wave resonance. ; Inlet pool swaying: fr = 0.1-0.5 Hz; Pump blade passing frequency harmonics.

[0399] B4. Frequency zone setting is prohibited:

[0400] For each resonant frequency fres,k, establish a forbidden interval:

[0401]

[0402] Among them: security boundary Select 1.5-2.5Hz (open channel system has low damping, and the restricted area width is larger than the pipe).

[0403] In the optimization constraints of step 4, the following should be added: the pump frequency fi must not fall into Fforbidden; the pump frequency difference |fi - fj| must not fall into Fforbidden (to avoid difference frequency resonance); the 2nd and 3rd harmonics of the pump frequency must not fall into Fforbidden.

[0404] The target flow rate Qtarget,i(fi, h) of the pump is calculated based on the following pump characteristic model and open channel hydraulic coupling:

[0405] Pump characteristic curve (corrected by similarity law):

[0406]

[0407]

[0408] Where: the efficiency ηi is fitted using a quadratic function:

[0409]

[0410] The optimal efficiency point flow rate (QBEP) varies linearly with frequency.

[0411] QBEP(fi) = QBEP,rated·(fi / frated)

[0412] Open channel hydraulic coupling:

[0413] The pump head needs to overcome: channel water level h, and pipeline friction loss. ,

[0414] Local loss .

[0415] Overall head balance equation:

[0416] Hi = h + hf + hm + hvp

[0417] Where: hvp is the difference between the pump installation elevation and the canal bottom elevation.

[0418] The target flow Qtarget,i is the solution to the following system of equations:

[0419]

[0420]

[0421] This set of equations demonstrates a strong coupling between pump characteristics and channel water level, and requires the use of Newton's iterative method to solve.

[0422] The calculation of unsteady flow in open channels employs the following two methods:

[0423] Method 1: Simplified calculation based on the method of characteristics:

[0424] For rectangular channels, a single-wave approximation is used:

[0425]

[0426] in: Let be the surface wave propagation speed.

[0427] Propagation time of water level disturbance caused by flow velocity regulation:

[0428]

[0429] .

[0430] Method 2: Detailed simulation based on the Saint-Venant equation:

[0431] Discretize the Saint-Venant equations using the Preissmann four-point implicit difference scheme:

[0432]

[0433]

[0434] :

[0435] .

[0436] Boundary conditions: Upstream: Downstream: The hQ relationship is determined by the gate equation.

[0437] Calculation outputs: water level time history h(x,t) for each cross-section, flow velocity time history v(x,t) for each cross-section, Froude number time history Fr(x,t), and flow velocity nonuniformity time history for the measurement section. .

[0438] Risk assessment: If max[Fr(x,t)] > 0.85, there is a risk of hydraulic jump. If the water level changes too quickly, The flow field stability is insufficient.

[0439] Method 2 is preferred, while Method 1 can be used for simple rectangular channels where the flow rate adjustment range is less than 20%.

[0440] The optimization problem can be solved using the following three algorithms:

[0441] Algorithm 1: Constrained Particle Swarm Optimization (PSO):

[0442] Particle definition: Each particle represents a frequency combination f = [f1, f2, ..., fn].

[0443] Speed ​​updates:

[0444]

[0445] Location update:

[0446]

[0447] Constraint handling: The penalty function method is used, and the augmented objective function is defined as follows:

[0448]

[0449] Where M is the large penalty factor, with a value of 10^6.

[0450] Parameter settings: Number of particles: 30-50 (open channel systems have more constraints and require more particles); Number of iterations: 40-60 generations; Inertia weight: w = 0.9→0.4 (linearly decreasing); Learning factor: c1 = c2 = 2.0.

[0451] Algorithm 2: Multi-objective genetic algorithm (NSGA-II):

[0452] When simultaneous optimization of energy consumption, accuracy, and stability is required, NSGA-II is adopted: Encoding: Integer encoding, gene = frequency level index; Fitness: Multi-objective vector. Pareto sort: Non-dominated sort + crowding distance.

[0453] Genetic operators: Selection: Tournament selection; Crossover: Simulated binary crossover (SBX), probability 0.9; Mutation: Polynomial mutation, probability 0.1; Output: Multiple solutions on the Pareto front, selected by the operator according to the verification requirements.

[0454] Algorithm 3: Sequential Quadratic Programming (SQP):

[0455] Suitable for small-scale problems with ≤3 pumps, with fast convergence speed;

[0456] Solve the quadratic programming subproblem in each iteration:

[0457]

[0458]

[0459] (Linearization constraints);

[0460] Update the Hessian matrix approximation using the BFGS formula.

[0461] Optimal algorithms: NSGA-II for precision mode, PSO for energy-saving mode, and SQP for real-time response.

[0462] The method further includes the following adaptive parameter update steps:

[0463] Periodic parameter calibration (performed every 48-168 hours of operation):

[0464] Collect operational log data, including frequency commands, measured flow rate, measured water level, and measured flow velocity distribution.

[0465] Update dynamic model parameters using system identification algorithms Extended Kalman Filter (EKF) or Unscented Kalman Filter (UKF) is used to fit the measured water level response curve with the model prediction residual.

[0466] Update the Manning roughness coefficient n: extract (Qsteady, hsteady, vsteady) from steady-state data; calculate it using the Manning formula: ; Detect n-value drift (such as n increasing due to algal attachment).

[0467] Update pump characteristic curve parameters H0, rated, KH, rated, ηmax: extract steady-state operating points (Qi, Hi, fi) from operating data; fit characteristic curve equations using multivariate nonlinear regression; detect impeller wear (alarm when H0 decreases by >5%).

[0468] Update the stability map M: Recalculate the ASI values ​​of key areas (commonly used calibration flow points) using new parameters; adopt an incremental update strategy, only recalculating areas with changes exceeding 15%; retain the original map as a backup and gradually replace it.

[0469] Online learning optimization (optional): Record each verification task (Qd, n, f) and verification results (uncertainty, repeatability); establish a decision-outcome database.

[0470] Stepwise optimization using reinforcement learning (Q-learning or Actor-Critic): State space: ( Action space: (n, f1, f2, ..., fn); Reward function: R = -uncertainty- 0.1·energy consumption + 10·(ASI>0.15).

[0471] Online correction of channel hydraulic characteristics:

[0472] Detecting siltation at the bottom of the canal: If the water level rises systematically under the same flow rate, it is inferred that the cross-sectional area of ​​the water passage has decreased;

[0473] Corrected effective cross-sectional area Recalculate the Froude number threshold.

[0474] Detect changes in gate characteristics: Identify the gate flow coefficient Cd and fit it to the measured Qout-h relationship; detect gate corrosion or deformation (alarm when Cd decreases by >10%).

[0475] Application examples:

[0476] (I) Application of the method described in this embodiment in the calibration system of open channel flowmeters:

[0477] The system configuration includes: Pump set: 3-6 centrifugal pumps in parallel, with a rated flow rate of 50-500 L / s and a head of 10-50 m per pump. Water storage tank: 50-200 m³ capacity, ensuring continuous calibration for more than 30 minutes.

[0478] Open channel structure: Flow straightening section: 5-10m in length, equipped with guide grids or flow straightening screens; Flow stabilization section: 15-30m in length; Measurement section: 10-20m in length, meeting the straight section requirements specified in ISO 4359; Outflow section: equipped with electric gates or overflow weirs.

[0479] The installation location of the flow meter under test is: the middle of the measuring section, upstream ≥ 10 times the channel width, downstream ≥ 5 times the channel width.

[0480] Control objective: Flow velocity uniformity in the measurement section. (Meets the requirements of JJG 164-2000 Class I equipment); Water level stability: ; Traffic stability: Froude number range: 0.3 < Fr < 0.6 (avoid critical flow and rapid flow); verification uncertainty: U ≤ 0.5% (k=2).

[0481] Verification process:

[0482] 1. Based on the flow rate being tested, query the optimal number of pumps n* and frequency f* from the stability map;

[0483] 2. Perform pump start-up and frequency regulation, and synchronously regulate the downstream gate;

[0484] 3. Wait for stabilization (ASI > 0.15 and fluctuation < threshold for 30 consecutive seconds);

[0485] 4. Perform verification using the standard table method or the volumetric method;

[0486] 5. Record calibration data and system operating parameters. );

[0487] 6. Automatically generate calibration certificates.

[0488] Special requirements:

[0489] The stability map M needs to cover the entire range of the flow meter under test (Qmin to Qmax);

[0490] The resonant frequency identification is performed once per quarter, during the intervals between calibration tasks.

[0491] Frequency adjustment is prohibited during the verification process; maintain constant operating conditions.

[0492] In low-temperature environments (water temperature <5°C), the Manning coefficient n and wave velocity c need to be corrected.

[0493] (II) Application of the method in irrigation canal flow monitoring:

[0494] The system configuration for irrigation canal flow monitoring includes:

[0495] Pump set: 2-4 mixed flow pumps or axial flow pumps, single pump flow rate 200-2000 L / s; Channel type: trapezoidal or rectangular cross section, bottom width 2-10m, design water depth 0.5-2m; Monitoring equipment: ultrasonic water level gauge + ADCP velocity profiler.

[0496] Control objective:

[0497] Irrigation flow rate accuracy: ;

[0498] To avoid channel overflow: h < 0.9 hdesign;

[0499] To prevent erosion of the canal bottom: v < vcritical (determined according to soil type, generally <1.5 m / s);

[0500] Energy saving priority: w1 = 0.7, w2 = 0.2, w3 = 0.1.

[0501] The ASI threshold The value is set to 0.12-0.15, lower than the calibration device (0.15-0.18). Because irrigation allows for greater flow fluctuations, the Froude number constraint is relaxed to 0.2 < Fr < 0.8.

[0502] Increase monitoring to prevent clogging:

[0503] When local blockage (abnormal flow velocity distribution) is detected, the flow rate is increased to flush the blockage, and ASI optimization is paused.

[0504] Seasonal parameter adjustments: High water season (high flow demand): Prioritize multiple pumps to reduce Fr; Low water season (low flow demand): Prioritize fewer pumps to improve single pump efficiency; Peak irrigation season: Force ASI > 0.15 to avoid water supply interruptions.

[0505] (III) Application of the method in hydrological station flow measurement:

[0506] The system configuration includes: pump set: 1-3 submersible pumps for artificial flow control (natural river channel + pump station water supply mode); test section: natural river channel or artificial flume, meeting the GB 50179-2015 standard; test method: velocity area method or buoy method.

[0507] Control objective: Maintain stable flow rate at the test section. Duration ≥ 10 minutes; water level fluctuation: Froude number: 0.2 < Fr < 0.5 (natural river channels are mostly slow-flowing).

[0508] Special constraints: Pump start-up and shutdown must be avoided during the testing period (to prevent disturbance); frequency adjustment should be prioritized to avoid start-up and shutdown; the ASI threshold can be reduced to 0.10 during low flow periods at night.

[0509] Coupled with natural inflow: Real-time monitoring of upstream inflow Qnatural; pump makeup flow rate Qpump = Qtarget – Qnatural; when Qnatural fluctuation > 10%, rapid frequency adjustment is triggered.

[0510] Safety measures during the flood season: When the water level exceeds the warning value, all pumps shall be shut down; when the Froude number approaches 1.0 (critical flow), adjustment shall be prohibited and the status quo shall be maintained; when a rainstorm warning is issued, the water level shall be lowered to a safe range in advance.

[0511] Application examples:

[0512] The invention will now be described in further detail with reference to theoretical analysis. All data in this specification are based on theoretical simulations.

[0513] Theoretical simulation case: Rectangular open channel flow rate calibration device.

[0514] I. System Parameter Settings:

[0515] Open channel geometric parameters: Channel type: rectangular cross section, bottom width B = 1.5 m, design water depth hdesign = 0.8 m, channel bottom slope S = 0.001 (1‰), Manning roughness coefficient n = 0.013 (smooth concrete), straightening section length 8 m, steady flow section length 20 m, measurement section length 5 m.

[0516] Pump set parameters (based on standard pump sample): Pump model IS80-65-160 single-stage centrifugal pump, quantity of pumps 4 (3 in operation, 1 standby), rated flow rate Qrated = 200 L / s = 0.2 m³ / s, rated head Hrated = 32 m, rated speed nrated = 2900 rpm (frequency 50Hz), rated power Prated = 22 kW, rated efficiency ηrated = 72%.

[0517] Pump characteristic curve fitting parameters (based on sample curve): H0,rated = 35.8 m; ; = 74%; a = 2.5 (efficiency curve parameter); QBEP,rated = 0.22 m³ / s.

[0518] Pipeline parameters: Pipe diameter D = 0.2 m, pipe length L = 15 m, friction coefficient f = 0.02, local resistance coefficient .

[0519] Verification flow range: Qmin = 50 L / s (0.05 m³ / s); Qmax = 600 L / s (0.6 m³ / s); Coverage velocity range: v = 0.04-0.5 m / s.

[0520] II. Dynamic Modeling:

[0521] Based on the state-space dynamic equations of the embodiment, a 5-dimensional state-space model (with 3 pumps in operation) is established:

[0522] State variables: x = [Q1, Q2, Q3, h, v]T

[0523] Dynamic equations:

[0524] dQ1 / dt = (1 / 1.5)[Qtarget,1(f1,h) - Q1] - 0.08[(Q1-Q2)+(Q1-Q3)]

[0525] dQ2 / dt = (1 / 1.5)[Qtarget,2(f2,h) - Q2] - 0.08[(Q2-Q1)+(Q2-Q3)]

[0526] dQ3 / dt = (1 / 1.5)[Qtarget,3(f3,h) - Q3] - 0.08[(Q3-Q1)+(Q3-Q2)]

[0527] dh / dt = (1 / 1.2)[(Q1+Q2+Q3) - Qout(h)]

[0528] dv / dt = (1 / 20)[(Q1+Q2+Q3) / 1.2 - v] - (9.81 / 20)(h - 0.8) - 1.82•v²

[0529] Parameter descriptions: τp = 1.5s (inertial time constant estimated based on pump sample); kc = 0.08 (coupling coefficient of 3 pumps in parallel, greater than 2 pumps); Ac = B•h ≈ 1.2 m² (calculated based on an average water depth of 0.8m); the last term coefficient 1.82 = gn² / ( ), where R ≈ 0.4m.

[0530] Outflow calculation (rectangular thin-walled weir):

[0531] .

[0532] III. Stability Map Generation (Theoretical Calculation):

[0533] The calculation script was written using MATLAB. Parameter space definition: Number of pumps: Verification flow rate: Segmented sampling: 0.05-0.2 m³ / s: step size 0.02 m³ / s; 0.2-0.6 m³ / s: step size 0.05 m³ / s; Frequency combination: f {35, 37, 39, ..., 50} Hz.

[0534] For each combination (n, Qd, f):

[0535] Step 1: Calculate steady-state water level for uniform flow in an open channel:

[0536] Given: Qd, B, n, S;

[0537] Solution: ;

[0538] For rectangular channels: ;

[0539] Substituting the values, we obtain the nonlinear equation:

[0540]

[0541] The steady-state water level hss is obtained by solving the problem using Newton's iteration method.

[0542] Step 2: Calculate the steady-state velocity and Froude number:

[0543]

[0544]

[0545] If Fr < 0.2 or Fr > 0.8, mark the condition as infeasible and skip it.

[0546] Step 3: Calculate the operating point of each pump:

[0547] For each pump i, solve the coupling equations:

[0548]

[0549]

[0550] The last term, 1.0m, represents the pump installation elevation (above the bottom of the channel). The iterative solution yields Qi,ss.

[0551] Step 4: Construct the Jacobian matrix (5×5):

[0552] At the steady-state point xss = [Q1,ss, Q2,ss, Q3,ss, hss, vss];

[0553] Calculate each element using numerical differentiation, with a step size. ;

[0554] Key coupling terms:

[0555] (Increased pump flow rate → rising water level);

[0556] (Rising water level → increased outflow);

[0557] (Water level rises → flow velocity decreases, due to gravity);

[0558] Complex terms, coupled through pump characteristic curves.

[0559] Step 5: Eigenvalue Analysis:

[0560] Calculate eigenvalues ​​using MATLAB's eig() function ;

[0561] calculate .

[0562] Step 6: Estimation of flow velocity uniformity:

[0563] Using the logarithmic law assumption for longitudinal velocity distribution:

[0564]

[0565] Where: v* is the frictional velocity, κ=0.41 is the Karman constant, and y0 is the roughness height;

[0566] Calculation of non-uniformity using the 5-point method of cross-section:

[0567]

[0568] like , marked as insufficient precision.

[0569] Example of calculation results (Qd = 0.3 m³ / s, 3 pumps):

[0570]

[0571] Total calculations: 17 (flow points) × 2 (number of pumps) × 8³ (frequency combinations) ≈ 17,408 operating points;

[0572] Computation time: Approximately 4-6 hours (for a standard PC);

[0573] Generated map size: approximately 3MB.

[0574] IV. Resonance Frequency Identification (Theoretical Analysis):

[0575] The main resonance sources in open channel systems:

[0576] (1) Channel standing wave resonance:

[0577] Baseband: ;

[0578] Second harmonic: .

[0579] Note: These frequencies are very low and generally do not coincide with the pump frequency, but the difference frequency may be close.

[0580] (2) Sloshing of the inlet pool:

[0581] Baseband: ;

[0582] Where: d is the water depth, and W is the pool width.

[0583] Assume d = 2m, W = 5m:

[0584] fs ≈ 0.35 Hz.

[0585] (3) Pump blade passing frequency:

[0586] For a 6-blade impeller, at 50Hz:

[0587] fb = 50 × 6 / 60 = 5 Hz

[0588] However, this is a mechanical vibration frequency, which has little impact on the surface of the open channel.

[0589] Theoretically predicted forbidden frequency region:

[0590] [0, 0.5] Hz: Covers standing waves and inlet pool sloshing;

[0591] [4.5, 5.5] Hz: Blade passing frequency (the impact is small if a flexible connecting tube is used);

[0592] Differential frequency forbidden zone: In practical applications, frequency sweep testing is required for verification.

[0593] V. Simulation of Saint-Venant Equations (Theoretical Case):

[0594] Unsteady flow simulations were performed using HEC-RAS software (an open-source open channel hydraulic calculation software developed by the U.S. Army Corps of Engineers).

[0595] Simulation scenario: The flow rate increases from 0.3 m³ / s to 0.45 m³ / s (an increase of 50%).

[0596] Traditional method (quick adjustment, completed in 5 seconds):

[0597] Simulation settings: Spatial step size Time step (Satisfies CFL conditions); Initial conditions: , Fr0 = 0.48; Boundary conditions: Q(0,t) = 0.3 + 0.03•t (0≤t≤5s).

[0598] Simulation results:

[0599] t=2s: The upstream water level rises to h=0.78m, Fr=0.52 (still safe);

[0600] t=3s: Fluctuations occur in the middle reaches, with an amplitude of ±3cm;

[0601] t=4s: The downstream water level begins to rise, Fr=0.68 (approaching the critical point);

[0602] t=5s: Fr in the middle of the measurement section reaches 0.82 (exceeding the safe range);

[0603] t=6s: Local occurrence of Fr>0.9, a precursor to a hydraulic jump;

[0604] Peak flow rate nonuniformity This exceeds the verification requirements.

[0605] The method described in this embodiment (based on S5 hydraulic jump pre-verification, adjusted to 10 seconds):

[0606] Simulation settings: Boundary conditions: Q(0,t) = 0.3 + 0.015•t (0≤t≤10s; Simultaneously adjust the gate opening to maintain a slow rise in water level.

[0607] Simulation results:

[0608] t=5s: Water level h=0.74m, Fr=0.50 (smooth rise);

[0609] t=10s: Water level h=0.76m, Fr=0.54 (still within safe range):

[0610] The maximum value of Fr throughout the entire process is 0.58, far from the critical flow;

[0611] Peak flow rate nonuniformity It meets the verification requirements.

[0612] Key improvement: By predicting the Fr trajectory, the adjustment time was extended from 5 seconds to 10 seconds, thus avoiding the risk of a water jump.

[0613] VI. Multi-objective optimization simulation (PSO algorithm)

[0614] Optimization issues:

[0615] Qd = 0.4 m³ / s;

[0616] 3 pumps; min F = [F1, F2, F3]

[0617] f;

[0618] F1 = ΣPi (total power, kW);

[0619] (Flow rate deviation, m³ / s);

[0620] (Flow rate non-uniformity).

[0621]

[0622] Comparison of optimization results:

[0623]

[0624] Analysis: The traditional optimization scheme has the lowest energy consumption, but its ASI value of 0.11 indicates instability. =3.8% accuracy insufficient; the solution of this invention improves ASI by 82% with only a 3% increase in energy consumption. Reduced by 50%.

[0625] VII. Stability Monitoring and Emergency Control Simulation:

[0626] Simulated scenario: Disturbances occur during operation (such as gate jamming);

[0627] Initial state: f=[49,47,44]Hz, h=0.74m, v=0.36m / s, ASI=0.20;

[0628] t=0s: Gate jams, Qout drops suddenly by 20%.

[0629] t=2s: Water level begins to rise, h=0.78m

[0630] t=4s: Detected (Positive value, signs of instability)

[0631] t=4.1s: Emergency control is triggered: frequency adjustment is frozen; the historical buffer is queried and the safe state f=[48,46,43] at t=-5s is found; Fr=0.61 is detected (still safe), and the frequency is adjusted instead of the pump is stopped;

[0632] t=4.2s: The command is issued, and the frequency of the three pumps is reduced synchronously, with a ramp time of 15s;

[0633] t=19s: The frequency drops to [48,46,43], and h drops to 0.76m;

[0634] t=20s: ASI recovers to 0.17. (Stablize);

[0635] t=25s: The gate fault is manually resolved, and the system returns to normal;

[0636] In contrast, without emergency control, Fr will exceed 0.85 at t=10s, potentially leading to a hydraulic jump.

[0637] VIII. Parameter Adaptive Simulation:

[0638] Scenario: Siltation at the bottom of the canal caused the Manning coefficient to increase from n=0.013 to n=0.016.

[0639] Initial (no sedimentation):

[0640] When Qd = 0.3 m³ / s, the steady-state h = 0.72m; the model predicts h = 0.72m, with an error of 0%.

[0641] After 3 months of operation (accumulation):

[0642] When Qd = 0.3 m³ / s, the measured h = 0.78 m; the model predicts h = 0.72 m, with an error of -7.7%.

[0643] Execution parameter identification:

[0644] Collect 100 steady-state operating points (Qi, hi, vi);

[0645] Calculate using Manning's formula:

[0646]

[0647] Statistical mean: n_ident = 0.0158 ≈ 0.016;

[0648] Update the model parameter n = 0.016.

[0649] After the update:

[0650] The model predicts h = 0.77m with an error of -1.3%, and the accuracy has been restored.

[0651] IX. Verification Process Simulation

[0652] Flow meter under test: Electromagnetic flow meter, range 0-1000 L / s

[0653] Verification points: 100, 200, 400, 600 L / s (0.1, 0.2, 0.4, 0.6 m³ / s)

[0654] Verification point 1: Qd = 0.1 m³ / s

[0655] Step 1: Query the stability map:

[0656] Recommendation: n=2, f=[46,44]Hz, expected ASI=0.21.

[0657] Step 2: Start the pump:

[0658] Start pumps 1 and 2, and ramp up the frequency from 0 to the target in 10-second increments.

[0659] Step 3: Adjust the gate:

[0660] Based on the measured water level, the closed-loop regulating gate is adjusted to h=0.65m.

[0661] Step 4: Stability assessment:

[0662] .

[0663] Step 5: Standard table method verification:

[0664] Simultaneously collect readings from the standard flow meter and the flow meter under test for 60 seconds;

[0665] Verification point 4: Qd = 0.6 m³ / s.

[0666] Step 1: Query the stability map:

[0667] Recommendation: n=3, f=[50,49,48]Hz, expected ASI=0.17.

[0668] Step 2: Switch from 2 pumps to 3 pumps:

[0669] Start the pump at 3 to 40 Hz, and then increase the frequency in tandem after it stabilizes.

[0670] Step 3: Adjust the gate:

[0671] Fully open the gate and maintain h=0.82m.

[0672] Step 4: Stability assessment:

[0673] Fr=0.64 (close to the upper limit), but ASI=0.18 (safe).

[0674] Permitted verification:

[0675] Full-range calibration time: 45-60 minutes for traditional methods, 30-40 minutes for the method of this invention (due to its stability and faster calibration).

[0676] 10. Theoretical Comparison with Traditional Methods:

[0677]

[0678] Finally, it should be noted that the above is only used to illustrate the technical solution of the present invention and not to limit it. Although the present invention has been described in detail with reference to the preferred arrangement, those skilled in the art should understand that modifications or equivalent substitutions can be made to the technical solution of the present invention (such as the form of the detection channel, the application of various formulas, the order of steps, etc.) without departing from the spirit and scope of the technical solution of the present invention.

Claims

1. An open channel flow measurement system based on attractor stability, the open channel flow measurement system comprising: A water storage tank, a pump set with at least two pumps connected in parallel, a flow straightening section equipped with a guide grid or flow straightening grid, a flow stabilization section, a measurement section, and an outlet section equipped with an electric gate or overflow weir are sequentially arranged in an open channel; the pump set is connected to a detection controller, which is equipped with a flow measurement unit, a water level measurement unit, a flow velocity measurement unit, a wave monitoring unit, an edge computing unit, a communication network module, a clock synchronization unit, and an execution control module; the detection controller is connected to a cloud server, characterized in that... The flow measurement unit is installed in the outlet pipe of each pump; The water level measurement unit is installed at one location each upstream, middle, and downstream of the measurement section; The flow velocity measurement unit is installed at points 5-9 on the cross section of the measurement section; The edge computing unit includes: A Fakubi matrix calculation module using numerical differentiation or automatic differentiation is employed. An eigenvalue solving module using the QR decomposition algorithm or the Jacobi iterative algorithm is employed. A multi-objective optimization solution module employing genetic algorithm, particle swarm optimization algorithm, or NSGA-II algorithm is used. A simulation module for unsteady open channel flow based on the finite difference method or finite volume method of the Saint-Venant equations; A stability map storage module that stores pre-computed M(n, Qd) data in the form of a database or memory-mapped file; The cloud server provides historical data storage services for storing ASI time series, Froude number time series, flow rate distribution, and calibration results; it also provides stability map update services, optimizes M parameters based on long-term operating data, and automatically generates calibration reports. The remote monitoring visualization interface displays the channel longitudinal profile water level line, flow velocity vector field, and attractor trajectory.

2. The system according to claim 1, characterized in that... The edge computing unit includes: The fast Lyapunov exponent estimation module employs the following algorithm: Extract the most recent M state vectors {x(t1), x(t2), ..., x(tM)} from the historical buffer, where M = 80-150; Calculate the distance between adjacent states: ； Calculate the normalized distance growth rate: ； This estimated value serves as a real-time stability monitoring indicator. ASI Fast Approximation Module: Calculate the Jacobian matrix at the current state x0; The Arnoldi iteration only calculates the 4-6 eigenvalues ​​with the largest modulus; Use the formula: ； Calculation time: 25-40ms; Module for rapid assessment of flow velocity uniformity; The cross-sectional velocity distribution {v1, v2, ..., vm} is measured using the 5-point method or the 9-point method. Calculate the flow velocity non-uniformity: ； Real-time determination of whether it meets the requirements .

3. A multi-pump coordinated method using the open channel flow rate determination system based on attractor stability as described in claim 1, characterized in that, The method includes the following steps: Step 1, Dynamic Modeling of the Open Channel-Pump System: Based on the measured characteristic curves of each pump and the hydraulic characteristics of the open channel, establish a set of state-space dynamic equations: ； Where: Qi is the instantaneous flow rate of the i-th pump; τp is the pump hydraulic time constant, in seconds; Qtarget,i is the target flow rate determined by the frequency fi and the channel water level h in the calibration section; kc is the pump coupling coefficient; Ac is the channel cross-sectional area; Qout is the outflow flow rate, determined by the downstream gate opening and water level; v is the average flow velocity in the measurement section; L is the length of the steady flow section; g is the gravitational acceleration; h0 is the channel bottom elevation; n is the Manning roughness coefficient; R is the hydraulic radius; N is the number of currently operating pumps. Step 2, real-time evaluation of attractor stability: Calculate the Jacobian matrix at the current running state x0 = [Q1, Q2, ..., QN, h, v]T: ； Each element is calculated using numerical differentiation: ； The step size δ is taken as 10^(-6) times the value of the corresponding variable; Calculate all eigenvalues ​​{λ1, λ2, ..., λn} of the Jacobian matrix; Define the attractor stability index ASI: ； Where β is the weighting coefficient, with a value ranging from 0.05 to 0.15, and an optimal value of 0.10; ε is the regularization constant, with a value of 0.01; Step 3, Dynamic optimization of pump combination number: Query the pre-built stability map M(n, Qd), which stores the ASI values ​​for different combinations of pump number n and calibration flow rate Qd; Based on the current calibration flow rate requirement Qd, select the number of pumps that maximizes ASI. ； Pump start / stop switching is performed when any of the following conditions are met: The current ASI is < 0.12; The expected energy savings after the switch will exceed 8%; |n* - ncurrent| ≥ 1 and the switching has been running stably for more than 300 seconds; Froude number Exceeding the safe range [0.3, 0.7]; Step 4, Frequency Allocation Multi-Objective Optimization: After determining the number of pumps n*, solve the optimization problem: ； ； Step 5, Pre-verification of hydraulic transient risks: Input the current state and the optimized target frequency f* into the open channel unsteady flow calculation model; Solve using the Saint-Venant equations: ； ； Calculate the water level change amplitude Δhmax and the flow velocity change rate dv / dt during the flow regulation process: Assessing the risk of a hydraulic jump: like If Fr is predicted to be > 0.8, indicating a risk of hydraulic jump, the following measures should be taken: Extend the adjustment time to ; Prioritize adjusting the downstream gate opening, and gradually adjust it in conjunction with the pump frequency; Phased adjustment: First adjust one pump, and adjust the other pumps after it stabilizes; Determining wave propagation: Calculate the surface wave propagation speed If the flow rate adjustment rate dv / dt > c / L, it may cause wave interference, and the adjustment should be delayed. Step 6, Control command execution and feedback monitoring: Send frequency command fi* and corresponding ramp time Δt to each inverter: Synchronous adjustment of downstream gate opening To maintain water level stability; Real-time acquisition of data from flow rate, water level, and flow velocity sensors, with a sampling period of no more than 1 second; Add the newly collected data to the historical buffer, with the buffer length set to 50-200 sampling points; Calculate the short-term Lyapunov exponent estimate within the sliding window: ； Where: d(t) is the distance between adjacent trajectory points in phase space, and Twindow is the window duration, which is 5-10 times. ; When continuously detected Emergency stabilization control is triggered if the duration exceeds 5 seconds, or if the following anomalies are detected: Water level fluctuations exceeding ±5cm and frequency >0.5Hz; The Froude number has exceeded 0.75; Flow velocity non-uniformity in the measurement section >5%; Emergency control measures: Freeze all frequency adjustment commands immediately: Query the historical buffer to find the most recent state fsafe that satisfies ASI > 0.18 and Fr∈[0.4,0.6]. Asymptotically recover to fsafe with maximum ramp time; Synchronously adjust the downstream gate to the corresponding opening degree; Send an alarm signal to the monitoring system; Return to step 2 and continue executing the loop at 1-second intervals until one of the following stopping conditions is met, at which point the loop stops: The current testing sites are consistently meeting the standards. The duration of flow deviation <1% exceeds the calibration sampling time; The operator issues a stop command through the monitoring interface or the local control cabinet; Single-point stability timeout (failed to meet standard after more than 15 minutes); An emergency fault was detected (pump failure, water level exceeding limit, communication interruption).

4. The method according to claim 3, characterized in that... The stability map M(n, Qd) is generated through the following offline calculation steps: A1, Parameter space discretization: Define the computational grid: The number of pumps ranges from n to {1, 2, ..., Ntotal}. Verification flow range: Qd ∈ [Qmin, Qmax], sampling interval segmented according to the range: Low flow rate range (0-50 L / s): 5 L / s interval; Medium flow rate range (50-500 L / s): 20 L / s interval; High flow rate range (500-2000 L / s): 50 L / s interval; For each number of pumps n, generate a frequency combination set Fn, requiring f1 ≥ f2 ≥ ... ≥ fn, with a frequency step size of 2Hz; A2. Steady-state hydraulic calculation of open channels: For each parameter combination (n, Qd, f): Calculate steady-state water level using the formula for uniform flow in an open channel: ； Where: S is the slope of the canal bottom; Obtain the operating points (Qi, Hi) of each pump based on the pump characteristic curve; The velocity distribution in the measurement section is calculated using either the logarithmic law or the power law. Verify whether the Froude number Fr is within a safe range; A3. Linearization of dynamics: Construct the Jacobian matrix at the steady-state solution xss = [Q1,ss, Q2,ss, ..., hss, vss]T: ； Using the central difference scheme: ； Where: ej is the j-th standard basis vector, and the step size δ is the scale of the variables. times; For open channel systems, the key partial derivatives include: The water level-velocity coupling term reflects the effect of gravity. The downstream outflow characteristics depend on the gate type; : Pump operating point drift with water level changes A4. Stability Quantification Calculation: Calculate the eigenvalue spectrum {λ1, λ2, ..., λn} of the Jacobian matrix; Calculate the ASI value; Constraint violation conditions: If the water level exceeds [hmin, hmax], it is marked as infeasible; if the Froude number Fr < 0.2 or Fr > 0.8, it is marked as a risky condition; if the flow velocity non-uniformity in the measurement section is >3%, it is marked as low accuracy; if ASI < 0, it is marked as unstable; if the efficiency of any pump is <40%, it is marked as inefficient. A5. Data Storage and Indexing: The calculation results (n, Qd, f, ASI, Fr, Store the total (Ptotal, flags) in an SQLite database or an HDF5 file; Create a multidimensional index: Primary key: (n, Qd); Indexed fields: ASI, Fr, .

5. The method according to claim 4, characterized in that... The resonant frequency fres is obtained through the following online identification method: B1. Small-amplitude incentive injection: When the system is running stably (selecting a medium flow rate condition, Fr≈0.5), inject a sinusoidal frequency disturbance into the selected pump: fi(t) = fi,0 + A·sin(2πftest·t) Wherein: the reference frequency fi,0 is the current operating frequency, the amplitude A is 0.3-0.8Hz, the test frequency ftest is scanned from 0.05Hz to 5Hz, and the duration of each test frequency is 30-60 seconds; B2. System Response Measurement: Simultaneously acquire the following at a sampling rate of not less than 20Hz: water level signal hw(t) in the middle of the measurement section; average flow velocity signal v(t) in the measurement section; optional: water surface fluctuation signal η(t); For each test frequency, extract the steady-state oscillation data and calculate the frequency response function: Hh(ftest) = |FFT[hw(t)]|f=ftest / |FFT[fi(t)]|f=ftest Hv(ftest) = |FFT[v(t)]|f=ftest / |FFT[fi(t)]|f=ftest Record the amplitudes |Hh|, |Hv| and the phases ∠Hh, ∠Hv; B3. Resonance Peak Identification: Detect frequency points fr that meet one of the following conditions in the frequency sweep results: Water level resonance criteria: amplitude is a local maximum, and |Hh(fr)| > 3×mean(|Hh|); half-power bandwidth < 0.2Hz (open channel resonance Q value is low); phase jump > 45° near fr; Velocity resonance criterion: ; accompanied by increased water surface ripples; The identified fr is recorded as the system's natural frequency. Typical resonance sources include: channel standing wave resonance: fr = c / (2L), c=(gh)^0.5; inlet pool sloshing: fr = 0.1-0.5 Hz; pump blade passing frequency harmonics; B4. Frequency zone setting is disabled: For each resonant frequency fres,k, establish a forbidden interval: ； The safety boundary Δf is taken as 1.5-2.5Hz; The following constraints are added to the optimization constraints in step S4: the pump frequency fi must not fall into Fforbidden; the pump frequency difference |fi - fj| must not fall into Fforbidden; and the 2nd and 3rd harmonics of the pump frequency must not fall into Fforbidden.

6. The method according to claim 5, characterized in that... The target flow rate Qtarget,i(fi, h) is calculated based on the following pump characteristic model and open channel hydraulic coupling: Pump characteristic curve: ； ； The efficiency ηi is fitted using a quadratic function: ； The optimal efficiency point flow rate (QBEP) varies linearly with frequency. QBEP(fi) = QBEP,rated·(fi / frated); Open channel hydraulic coupling: ； Overall head balance equation: Hi = h + hf + hm + hvp Where: hvp is the difference between the pump installation elevation and the canal bottom elevation; The target flow Qtarget,i is the solution to the following system of equations: ； 。 7. The method according to claim 6, characterized in that... It also includes the following adaptive parameter update steps: Periodic parameter calibration: Collect operational log data, including frequency commands, measured flow rate, measured water level, and measured flow velocity distribution; Update dynamic model parameters using system identification algorithms : Use either an extended Kalman filter (EKF) or an unscented Kalman filter (UKF); Fit the measured water level response curve to the model prediction residual; Update Manning roughness coefficient n: Extract (Qsteady, hsteady, vsteady) from steady-state operating data; Using Manning's formula to calculate inversely: ; Detect n-value drift; Update the pump characteristic curve parameters H0,rated, KH,rated, ηmax: Extract steady-state operating points (Qi, Hi, fi) from the operating data; The characteristic curve equation was fitted using multivariate nonlinear regression. Inspect impeller wear; Update stability map M: Recalculate the ASI values ​​for the critical region using the new parameters; An incremental update strategy is adopted, and only areas where changes exceed 15% are recalculated. The original map is retained as a backup, and the replacement is performed gradually. Online learning optimization: Record each verification task (Qd, n, f) and the verification result: Establish a decision-outcome database; Gradual optimization using reinforcement learning: State space: (Qd, ASI, Fr, σv); Action space: (n, f1, f2, ..., fn); Reward function: R = -uncertainty - 0.1·energy consumption + 10·(ASI>0.15); Online correction of channel hydraulic characteristics: Detecting siltation at the bottom of the canal: If the water level rises systematically under the same flow rate, it is inferred that the cross-sectional area of ​​the water passage has decreased. Corrected effective cross-sectional area ; Recalculate the Froude number threshold; Detect changes in gate characteristics: Identify the gate flow coefficient Cd and fit it to the measured Qout-h relationship; detect gate corrosion or deformation.

8. The method according to claim 7, characterized in that... The following method is used to calculate the unsteady flow in the open channel: Simplified calculation based on the method of characteristics: For rectangular channels, a single-wave approximation is used: ； in: The surface wave propagation speed; Propagation time of water level disturbance caused by flow velocity regulation: ； 。 9. The method according to claim 7, characterized in that... The following method is used to calculate the unsteady flow in the open channel: Detailed simulation based on Saint-Venant's equations: Discretize the Saint-Venant equations using the Preissmann four-point implicit difference scheme: ； ； ： ； Boundary conditions: Upstream: Q(0,t) = ΣQi(t); Downstream: The hQ relationship is determined by the gate equation; Calculation output: water level time history h(x,t) for each cross-section; flow velocity time history v(x,t) for each cross-section; Froude number time history Fr(x,t); flow velocity nonuniformity time history σv(t) for the measurement section; Risk assessment: If There is a risk of hydraulic jump; if The water level changes too quickly; if The flow field stability is insufficient.

10. The method according to claim 8 or 9, characterized in that... The optimization problem is solved using the following algorithm: Constrained particle swarm optimization; Particle definition: Each particle represents a frequency combination f = [f1, f2, ..., fn]; Speed ​​updates: ； Location update: ； Constraint handling: The penalty function method is used, and the augmented objective function is defined as follows: ； Where: M is the large penalty factor, with a value of 10^6; Parameter settings; Number of particles: 30-50; Number of iterations: 40-60 generations; Inertia weight. ; Learning factors: .

11. The method according to claim 8 or 9, characterized in that... The optimization problem is solved using the following algorithm: Multi-objective genetic algorithm method: When it is necessary to optimize energy consumption, accuracy, and stability simultaneously, NSGA-II should be used. Encoding: Integer encoding, gene = frequency rank number; Fitness: Multi-objective vector ; Pareto sort: Non-dominated sort + crowding distance; Genetic Operators: Selection: Tournament Selection; Crossover: Simulated binary crossover, probability 0.9; Mutation: Polynomial mutation, probability 0.1; Output: Multiple solutions on the Pareto front.

12. The method according to claim 8 or 9, characterized in that... The optimization problem is solved using the following algorithm: Sequential quadratic programming method: Suitable for small-scale problems with ≤3 pumps, with fast convergence speed; Solve the quadratic programming subproblem in each iteration: ； Update the Hessian matrix approximation using the BFGS formula; Optimal algorithms: NSGA-II for precision mode, PSO for energy-saving mode, and SQP for real-time response.

Images