Method for river water and sediment calculation based on muskingum method

Abstract

The application discloses a river water and sediment calculation method based on the Muskingum method, and belongs to a general method in the water conservancy field. The method comprises the following steps: (1) target river section selection; (2) hydrological data collection; (3) flow and sediment transport rate data sequence compilation; (4) determining the value of the flow calculation coefficient; (5) using the Muskingum flow calculation formula to calculate the flow process; (7) using the Muskingum sediment transport rate calculation formula to calculate the sediment transport process. The application expands the Muskingum sediment transport rate calculation formula on the basis of the traditional Muskingum flow calculation formula, and obtains the Muskingum sediment transport rate calculation formula of the river section by simultaneously solving the river section sediment balance equation and the sediment volume storage equation. The application synchronously realizes the fast and effective calculation of the water and sediment processes, has strong applicability and generalization, and can be used as a fast and effective hydrological method for the water and sediment calculation of the multi-sediment river.

Description

Technical Field

[0001] This invention pertains to methods for calculating river sediment, and in particular relates to a method for calculating river water and sediment based on the Muskingen method. Background Technology

[0002] Suspended sediment transport calculation is an important research subject for sediment-rich rivers, exemplified by the Yellow River. Currently, rapid hydrological calculation methods are lacking for sediment calculation during flood evolution in large rivers like the Yellow River. Hydrodynamic-based numerical models of river sediment transport are the primary tools for calculating sediment transport in the Yellow River's main channel and tributaries; however, these models are complex, computationally intensive, and subject to numerous boundary value conditions, making them unsuitable for channel sediment transport calculations in large-basin hydrological models. Existing research has introduced artificial neural networks into river sediment transport calculations, reducing reliance on topography compared to physical process-based hydrodynamic models, but lacking information on sediment transport mechanisms. For some rivers lacking measured topographic data, when numerical models cannot be used to calculate sediment transport processes, finding accurate hydrological methods to calculate sediment transport becomes crucial. In addition, when applying hydrological models to calculate sediment yield and transport in the Yellow River Basin, hydrodynamic models or hydrological methods are generally used to calculate river network confluence. In gullies, balance or unbalance sediment transport calculations are mostly based on the assumption of steady flow, which is theoretically contradictory to the use of unsteady flow calculations for water flow. Summary of the Invention

[0003] Purpose of the invention: The purpose of this invention is to provide a method for calculating river water and sediment based on the Muskingen method.

[0004] Technical solution: The present invention provides a method for river water and sediment calculation based on the Muskingen method, comprising the following steps:

[0005] S1. Determine the target river section to be calculated, select appropriate hydrological stations as the upstream and downstream hydrological stations of the target river section, and determine the river length L between the two hydrological stations.

[0006] S2. Collect historical hydrological data from upstream and downstream hydrological stations, and compile information on flow and sediment transport rates for the same period at both stations, including the upstream hydrological station flow sequence qu, upstream flow time series tuq, upstream hydrological station sediment transport rate sequence qSu, upstream sediment transport rate time series tus, downstream hydrological station flow sequence qd, downstream flow time series tdq, downstream hydrological station sediment transport rate sequence qSd, and downstream sediment transport rate time series tud; determine the flood propagation time at the two hydrological stations, and select the time interval of the compiled data as an approximation of the flood propagation time or a divisor of its approximation;

[0007] S3. Use linear interpolation to standardize the flow and sediment transport rate data sequences of upstream and downstream hydrological stations, so that the data sequences of upstream and downstream hydrological stations are at the same time interval, at the same moment, and at the same time period.

[0008] S4. Use trial-and-error or least squares method to determine the values ​​of flow calculation coefficients C0, C1, and C2;

[0009] S5. The flow process is calculated using the Muskingum method flow calculation formula to obtain the flow calculation sequence of downstream hydrological stations;

[0010] S6. The sediment transport process was calculated using the Muskingan method sediment transport rate calculation formula to obtain the sediment transport rate calculation sequence for downstream hydrological stations; the Muskingan method sediment transport rate calculation formula is as follows:

[0011] Q Spd (t i+1 ) = C S,0 Q Su (t i+1 )+C S,1 Q Su (t i )+C S,2 Q Sd (t i )+C S,3 S S (t i ), i = 1, 2, ..., N-1

[0012] Among them, C S,0 ~C S,3 C is the coefficient for calculating the sediment transport rate. S,0 C S,1 C S,2 The values ​​of C0, C1, and C2 are the same as those of the flow calculation coefficients, respectively. S,3 =C S,0 +C S,1 S S (t i ) for t i Riverbed erosion and sedimentation source; Q Su (t i+1 ), Q Su (t i ) are respectively t i+1 t i Sediment transport rate at upstream hydrological stations at any given time; Q Sd (t i ) for t i Sediment transport rate at downstream hydrological stations at any given time; Q Spd (t i+1 ) for t i+1The calculated sediment transport rate values ​​at downstream hydrological stations at given time points constitute the sediment transport rate calculation sequence for downstream hydrological stations as follows: QSPD={Q Spd (t2), Q Spd (t3), ..., Q Spd (t N )};

[0013] The accuracy of the method is determined by comparing it with the standard sequence of sediment transport rates at downstream hydrological stations, according to the following formula for calculating the coefficient of determination:

[0014]

[0015] in, The coefficient of determination of sediment transport rate is the formula for sediment transport rate calculation using the Muskingen method. The closer the value of Q is to 1, the higher the calculation accuracy. Spd (t i ) for t i Calculated sediment transport rate at downstream hydrological stations at any given time; This is the average value of the sediment transport rate standard sequence QSD for downstream hydrological stations, excluding the first term.

[0016] Furthermore, step S2, which involves compiling the flow and sediment transport rate information of the two hydrological stations during the same period, specifically includes:

[0017] The upstream hydrological station flow sequence qu is represented as: qu={q u (t uq 1 ), q u (t uq 2 ), q u (t uq 3 ), ..., q u (t uq n1 )};

[0018] The corresponding upstream flow time series tuq is represented as: tuq={t uq 1 , t uq 2 , t uq 3 , ..., t uq n1};

[0019] Among them, t uq 1 t uq 2 t uq 3 t uq n1 These represent the times of the flow sequence at the upstream hydrological stations at the 1st, 2nd, 3rd, and nth consecutive moments, respectively. u (t uq 1 ), q u (t uq 2 ), q u (t uq 3 ), q u (tuq n1 ) represent t respectively uq 1 t uq 2 t uq 3 t uq n1 The flow rate at the upstream hydrological station at a given time, with a flow sequence length of n1;

[0020] The sediment transport rate sequence qSu of upstream hydrological stations is represented as: qSu={q Su (t us 1 ), q Su (t us 2 ), q Su (t us 3 ), ..., q Su (t us n2 )};

[0021] The corresponding upstream sediment transport rate time series tus is represented as: tus={t us 1 , t us 2 , t us 3 , ..., t us n2};

[0022] Among them, t us 1 t us 2 t us 3 t us n2 These represent the times of the sediment transport rate sequence at the upstream hydrological stations, specifically the times of the 1st, 2nd, 3rd, and n2nd points, respectively. Su (t us 1 ), q Su (t us 2 ), q Su (t us 3 ), q Su (t us n2 ) represent t respectively us 1 t us 2 t us 3 t us n2 The sediment transport rate of upstream hydrological stations at any given time, and the length of the sediment transport rate sequence of upstream hydrological stations is n².

[0023] The downstream hydrological station flow sequence qd is expressed as: qd={q d (t dq 1 ), q d (t dq 2 ), q d (t dq 3 ), ..., q d (t dq n3 )};

[0024] The corresponding downstream flow time series tdq is represented as: tdq={t dq 1 , t dq 2 , t dq 3 , ..., tdq n3};

[0025] Among them, t dq 1 t dq 2 t dq 3 t dq n3 q represents the times of the flow sequence at downstream hydrological stations at the 1st, 2nd, 3rd, and nth consecutive moments. d (t dq 1 ), q d (t dq 2 ), q d (t dq 3 ), q d (t dq n3 ) represent t respectively dq 1 t dq 2 t dq 3 t dq n3 The flow rate at downstream hydrological stations at given time, and the length of the flow rate sequence at downstream hydrological stations is n3;

[0026] The sediment transport rate sequence qSd at downstream hydrological stations is expressed as: qSd={q Sd (t ds 1 ), q Sd (t ds 2 ), q Sd (t ds 3 ), ..., q Sd (t ds n4 )};

[0027] The corresponding downstream sediment transport rate time series tud is expressed as: tud={t ud 1 , t ud 2 , t ud 3 , ..., t ud n4};

[0028] Among them, t ud 1 t ud 2 t ud 3 t ud n4 These represent the times of the sediment transport rate sequence at downstream hydrological stations, specifically the 1st, 2nd, 3rd, and nth times respectively. Sd (t ds 1 ), q Sd (t ds 2 ), q Sd (t ds 3 ), q Sd (t ds n2 ) represent t respectively ds 1 t ds 2 t ds 3 t ds n4 The sediment transport rate of downstream hydrological stations at a given time, and the length of the sediment transport rate sequence of downstream hydrological stations is n4;

[0029] The method for calculating the flood propagation time T is as follows:

[0030] Based on the compiled information on flow and sediment transport rate of the two hydrological stations during the same period, the corresponding times of the flood peaks at the upstream and downstream hydrological stations during the same flood event are identified. The time difference between the two is used to determine the flood propagation time of the event. The flood propagation time is different for different flood events. The average value is selected to determine the flood propagation time T between the two hydrological stations.

[0031] The time interval Δt for compiling the data is chosen as an approximation of the flood propagation time or a divisor of that approximation.

[0032] Furthermore, step S3 specifically includes:

[0033] Let T1 be the maximum value of the first term of the four time series of flow and sediment transport rate at upstream and downstream hydrological stations, T1 = max{t uq 1 , t us 1 , t dq 1 , t ds 1 Let T2 be the minimum value of the last term of the four time series, then T2 = min{t}. uq n1 , t us n2 , t dq n3 , t ds n4}, calculate (T2-T1) / Δt+1 and take its integer value as the length N of the unified time series; starting from t1=T1, determine the unified time series {t1, t2, t3, ..., t4} with the time interval Δt as the common difference. N}, where t1, t2, t3, t N Let t represent the times 1, 2, 3, and N in sequence. This time series is an arithmetic sequence, and its start and end points are all within the four time series ranges of flow and sediment transport rate at upstream and downstream hydrological stations. For any term t within the unified time series... i Let i take values ​​from 1 to N, and find t. i At the four time series positions of flow and sediment transport rate at upstream and downstream hydrological stations, if t i In the upstream flow time series tuq, t satisfies uq k <t i <t uq k+1 Where k is some integer, then Q u (t i )=q u (t uq k )+(q u (t uq k+1 )-q u (t uq k ))(t i -t uq k ) / (t uq k+1 -t uq k ), Q u(t i ) for t i The flow rate at upstream hydrological stations is used to obtain the upstream hydrological station flow rate specification sequence QU, which is represented as: QU={Q u (t1), Q u (t2), Q u (t3), ..., Q u (t N )}

[0034] Among them, Q u (t1), Q u (t2), Q u (t3), Q u (t N ) represent t1, t2, t3, t respectively. N Similarly, the sediment transport rate standard sequence QSU for upstream hydrological stations is obtained as follows:

[0035] QSU={Q Su (t1), Q Su (t2), Q Su (t3), ..., Q Su (t N )}

[0036] Among them, Q Su (t1), Q Su (t2), Q Su (t3), Q Su (t N ) represent t1, t2, t3, t respectively. N Sediment transport rate at upstream hydrological stations at all times;

[0037] Similarly, the downstream hydrological station flow specification sequence QD and sediment transport rate specification sequence QSD are respectively expressed as:

[0038] QD={Q d (t1), Q d (t2), Q d (t3), ..., Q d (t N )}

[0039] QSD={Q Sd (t1), Q Sd (t2), Q Sd (t3), ..., Q Sd (t N )}

[0040] Among them, Q d (t1), Q d (t2), Q d(t3), Q d (t N ) represent t1, t2, t3, t respectively. N Downstream hydrological station flow rate, Q Sd (t1), Q Sd (t2), Q Sd (t3), Q Sd (t N ) represent t1, t2, t3, t respectively. N Sediment transport rate at downstream hydrological stations at any given time.

[0041] Furthermore, in step S4, when using the trial-and-error method, the initial time t1 of the channel storage capacity WP(t1) is set to zero, and the values ​​of each term in the channel storage capacity sequence WP are calculated using the following formula: WP={WP(t1), WP(t2), WP(t3), ..., WP(t... N )},WP(t i ) = WP(t i-1 )+(Q u (t i )-Q d (t i ))Δt, i=2, 3,..., N;

[0042] Among them, WP(t1), WP(t2), WP(t3), WP(t N ) represent t1, t2, t3, t respectively. N The water storage capacity of the channel in the river section at any given time; WP(t) i ),WP(t) i-1 ) represent t respectively i t i-1 The water storage capacity of the river section at any given time; Q u (t i ), Q d (t i ) represent t respectively i The flow rates at upstream and downstream hydrological stations at specific times; the time interval Δt is the length of the calculation period.

[0043] The values ​​of each term in the reservoir discharge sequence QP of the river segment are calculated using the following formula for the discharge specification sequences QU and QD of the upstream and downstream hydrological stations: QP={QP(t1), QP(t2), QP(t3), …, QP(t… N )},QP(t i )=x Q u (t i )+(1-x)Q d (t i ), i = 1, 2, ..., N;

[0044] wherein, QP(t1), QP(t2), QP(t3), QP(t N ) respectively represent the indicating storage flow of the river reach at time t1, t2, t3, t N ; Q u (t i ), Q d (t i ) are respectively the flow of the upstream hydrological station and the flow of the downstream hydrological station at time t i , QP(t i ) is the indicating storage flow of the river reach at time t i ;

[0045] the flow specific gravity factor x is in the range of 0-1, and the set values are x1, x2, x3 respectively, which satisfy 0<x1<x2<x3<1. The indicating storage flow sequences QP1, QP2 and QP3 corresponding to x1, x2 and x3 are calculated by using the indicating storage flow calculation formula of the river reach. The relationships between the indicating storage flow sequences QP1, QP2, QP3 and the river reach channel storage water volume sequence WP are plotted respectively, and the x value that enables the closest linear relationship between the two is determined through observation and comparison. When the x value is x1, x4, x5 and x6 are taken, which satisfy 0<x4<x5=x1<x6<x2; when the x value is x2, x4, x5 and x6 are taken, which satisfy x1<x4<x5=x2<x6<x3; when the x value is x3, x4, x5 and x6 are taken, which satisfy x2<x4<x5=x3<x6<1. The indicating storage flow sequences QP4, QP5 and QP6 corresponding to x4, x5 and x6 are calculated by using the indicating storage flow calculation formula of the river reach. The relationships between the indicating storage flow sequences QP4, QP5, QP6 and the river reach channel storage water volume sequence WP are plotted respectively, and the x value that enables the closest linear relationship between the two is determined through observation and comparison. The above process is repeated in sequence until no further determination can be made through graphic observation and comparison; the slope value of the linear relationship between the indicating storage flow Q' and the river reach channel storage water volume W is the storage constant K. After K and x are determined through trial calculation, the flow routing coefficients C0, C1 and C2 need to be converted and determined by the following formula:

[0046]

[0047] further, when the least square method is used for calculation in step S4, the following matrix operation is directly used to solve and determine the values of flow routing coefficients C0, C1 and C2:

[0048]

[0049]

[0050] wherein, Y is an N-1 dimensional column vector composed of the standardized flow sequence QD of the downstream hydrological station, and X is an (N-1)×3 dimensional matrix composed of the standardized flow sequences QD and QU of the downstream hydrological station and the upstream hydrological station, X is a three-dimensional column vector composed of flow calculation coefficients. T Let X be the transrank matrix.

[0051] Furthermore, the Muskingan method flow calculation formula in step S5 is as follows:

[0052] Q pd (t i+1 )=C0Q u (t i+1 )+C1Q u (t i )+C2Q d (t i ), i = 1, 2, ..., N-1

[0053] Among them, Q u (t i+1 ), Q u (t i ) are respectively t i+1 t i Flow rate at upstream hydrological stations; Q d (t i ) for t i Downstream hydrological station flow rate; Q pd (t i+1 ) for t i+1 The calculated flow values ​​at downstream hydrological stations can form a flow calculation sequence for downstream hydrological stations: QPD={Q pd (t2), Q pd (t3), ..., Q pd (t N )};

[0054] The accuracy of the method is determined by comparing it with the flow specification sequence of downstream hydrological stations, according to the following formula for calculating the flow determination factor: in, Q is the flow determination factor calculated using the Muskingan method flow calculus formula. pd (t i ) for t i Calculated flow values ​​at downstream hydrological stations at any given time; The flow specification sequence QD for downstream hydrological stations is the average value excluding the first term; the flow determination coefficient. The closer the value is to 1, the higher the calculation accuracy.

[0055] Furthermore, in step S6, the source term is determined based on the relationship between the sediment transport capacity and sediment concentration of a specific river channel. When an unbalanced sediment transport method is adopted, we have:

[0056]

[0057] Wherein, S(t) i), S * (t i ) represent the average suspended sediment concentration and sediment-carrying capacity of the river section, respectively; ω represents the settling velocity of sediment particles; α represents the recovery saturation coefficient of suspended sediment; h(t) represents ... i ) represents the average water depth of the river cross section; L represents the river length between the two hydrological stations; A(t) represents the average water depth of the river cross section; L represents the average water depth of the river cross section between the two hydrological stations; A(t) represents the average water depth of the river cross section between the two hydrological stations. i () represents the cross-sectional area of ​​the river channel.

[0058] Alternatively, the source term can be calculated using empirical formulas based on the actual river section conditions. When the source term is not considered, S S (t i =0, and the formula for calculating sediment transport rate and flow rate is consistent.

[0059] The present invention provides a river sediment calculation system based on the Muskingan method, comprising:

[0060] The data acquisition and processing module is used to determine the target river section, select upstream and downstream hydrological stations for the target river section, and determine the river length between the two hydrological stations; collect historical hydrological data from the upstream and downstream hydrological stations, compile the flow and sediment concentration information of the two hydrological stations in the same period, determine the flood propagation time of the two stations, and select the time interval of the compiled data as an approximation or a divisor of the flood propagation time; compile the flow and sediment concentration (or sediment transport rate) data sequences of the upstream and downstream hydrological stations, and perform linear interpolation on the measured data to ensure that the data sequences of the upstream and downstream hydrological stations are at the same time interval, at the same moment, and in the same period.

[0061] The flow calculation coefficient calculation module is used to determine the value of the flow calculation coefficient using a trial-and-error method or the least squares method;

[0062] The flow process calculation module is used to calculate the flow process using the Muskingum method flow calculation formula to obtain the calculated value of the flow process at the downstream station.

[0063] The sediment transport process calculation module is used to calculate the sediment transport process using the Muskingen method sediment transport rate calculation formula to obtain the calculated value of the sediment transport process at the downstream station.

[0064] An apparatus of the present invention includes a memory and a processor, wherein:

[0065] Memory is used to store computer programs that can run on a processor;

[0066] A processor, used to execute, while running the computer program, the steps of a river sediment calculation method based on the Muskingen method described above.

[0067] The present invention provides a storage medium storing a computer program, which, when executed by at least one processor, implements the steps of the above-described method for river water and sediment calculation based on the Muskingen method.

[0068] Beneficial effects: Compared with the prior art, the significant technical effects of the present invention are as follows: The river sediment calculation method based on the Muskingen method can guarantee a certain accuracy and the calculation process is very simple. It can be applied to the sediment calculation of sediment-laden rivers with suspended sediment transport as the main feature. For hydrological forecasting using the Muskingen method for flood processes, this method can be used to predict the sediment transport process of floods. This model can be used to calculate the channel sediment transport process in large-scale watershed sediment models. Since the material transport equation represented by suspended sediment has universality, the established method is also applicable to the river material flux, including solutes, while avoiding the complex calculations of hydrodynamic models, thus having an advantage in computational efficiency. Attached Figure Description

[0069] Figure 1 This is a flowchart of the method of the present invention. Detailed Implementation

[0070] The present invention will now be described in detail with reference to the accompanying drawings and specific embodiments.

[0071] The principle of the method of this invention is as follows:

[0072] The Muskingen flow calculation algorithm obtains the flow calculation formula for a river segment by simultaneously solving the river segment water balance equation and the channel storage equation:

[0073] Water balance equation:

[0074] Storage equation: W=KQ′=K[xQ] u +(1-x)Q d ]

[0075] Muskingan method flow calculation formula:

[0076] Where W represents the channel storage capacity of the river section; Q u Q d These are the inlet and outlet flow rates of the river section, respectively; Q′ is called the storage capacity; K is the storage constant, which is equivalent to the propagation time of the flood wave between the two hydrological stations at the inlet and outlet of the river section; x is the flow weight factor, which is mainly related to the degree of flattening and deformation of the flood wave. These represent the inlet and outlet flows of the river segment at time step n+1, respectively. C0 and C1 represent the inlet and outlet flows of the river segment at time step n, respectively. Time steps n and n+1 correspond to times t and t+Δt, respectively, where Δt is the calculation period length. C0, C1, and C2 are flow calculation coefficients. The relationship between the flow calculation coefficients and the storage constant and flow weight factor is as follows:

[0077]

[0078] The sediment transport process of this invention mainly targets sediment transport where suspended matter is the primary mode of movement. It is assumed that the sediment transport rate within the river section changes linearly along the river, regardless of whether the transport rate increases or decreases. When the indicated storage flow rate and the channel storage capacity of the river section have a single relationship, an indicated storage sediment transport rate Q′ with the same weight as Q′ is adopted. S The amount of sediment stored in the channel of the river section, W S It can also be a single linear relationship, referred to as the Muskingen method sediment transport channel relationship or sediment transport channel equation:

[0079] W S =K S Q′ S =K S [x S Q S,u +(1-x S )Q S,d ]

[0080] Among them, W S K represents the amount of sediment stored in the channel of the river section. S Q′ is the sediment transport constant, equivalent to the propagation time of sediment between the two hydrological stations at the river's inlet and outlet. S This is called the sediment transport rate; x S Here, K represents the sediment density factor. It is assumed that the sediment transport capacity constant and sediment density factor are consistent with the storage constant and flow rate factor used in the river reach flow calculation, i.e., K. S =K, x=x S Q S,u Q S,d These represent the sediment transport rates at the river inlet and outlet, respectively. The sediment within the river section also satisfies mass balance; the change in sediment mass per unit time equals the difference between the sediment transport rates at the river inlet and outlet, plus the riverbed scouring and deposition source term. Therefore, the sediment balance equation is:

[0081]

[0082] Among them, S S For the source of river erosion and siltation;

[0083] Solving the sediment balance equation and sediment transport channel equation simultaneously yields the Muskingan method formula for calculating sediment transport rate:

[0084]

[0085] in, These represent the sediment transport rates at the river inlet and outlet of the river segment at time step n+1, respectively. These represent the sediment transport rates at the river inlet and outlet of the river segment at time step n; C represents the riverbed scouring and sedimentation source term at time step n; S,0 ~C S,3 C is the coefficient for calculating the sediment transport rate. S,0 C S,1 C S,2 The values ​​of C0, C1, and C2 are the same as those of the flow calculation coefficients, respectively. S,3 =C S,0 +C S,1 .

[0086] The calculation formula for flow rate and sediment transport rate established in this invention can be performed in segments, which is known as the segmented Muskingen method. The source term is determined based on the relationship between the sediment transport capacity and sediment concentration of a specific river channel. When using the unbalanced sediment transport method, we have:

[0087]

[0088] Among them, S, S * Here, ω represents the average suspended sediment concentration and sediment-carrying capacity of the river section; α is the settling velocity of sediment particles; h is the average water depth of the river cross-section; L is the river length; and A is the cross-sectional area of ​​the river. The source term can also be calculated using empirical formulas based on the actual river section conditions. When the source term is not considered, the formulas for sediment transport rate and flow rate are consistent. In this method, the water-sediment propagation time corresponds to the velocity of the water flow and sediment. This method uses the same propagation time as the water flow when calculating the sediment transport process, which means that this method assumes that the two propagate synchronously. When the propagation velocities of the two are inconsistent, the sediment peak lags behind the flood peak, and it is necessary to select the sediment propagation time for calculation; in this case, the sediment specific gravity factor x... S Still consistent with the flow proportion factor x S =x, the sediment propagation time is taken to be consistent with the time of water particle movement, that is, the speed of sediment moving downstream is equal to the speed of water flow, and the sediment propagation time is calculated using the following formula:

[0089] Ks=2L / (V um +V dm )

[0090] Among them, K S It refers to the time it takes for sediment to travel, V um V dm These are the average flow velocities at the peak flood times of the upstream and downstream hydrological stations, respectively. The following coefficients are used to calculate the sediment transport rate:

[0091]

[0092] C S,3 =C S,0 +C S,1

[0093] The one-dimensional suspended sediment transport equation and the kinematic wave equation have the same equation form and approximate wave velocities. The Muskingen flow calculation method, based on the kinematic wave equation, is also applicable to the suspended sediment transport equation. This invention's river sediment calculation method based on the Muskingen method ensures a certain level of accuracy and has a very simple calculation process, making it applicable to sediment calculations in sediment-rich rivers where suspended sediment transport is dominant. For hydrological forecasting using the Muskingen method for flood events, this method can predict flood sediment transport processes. This model can be used to calculate channel sediment transport processes in large-scale watershed sediment models. Because the suspended sediment transport equation has universality, the established method is also applicable to river material fluxes, including solutes, while avoiding the complex calculations of hydrodynamic models, thus offering advantages in computational efficiency.

[0094] like Figure 1 As shown, the present invention provides a river water and sediment calculation method based on the Muskingen method, comprising the following steps:

[0095] S1. Determine the target river section to be calculated, select appropriate hydrological stations as the upstream and downstream hydrological stations of the target river section, and determine the river length L between the two hydrological stations.

[0096] S2. Collect historical hydrological data from upstream and downstream hydrological stations, and compile information on flow and sediment transport rates for the same period at both stations. This includes the upstream station flow sequence qu, upstream flow time series tuq, upstream hydrological station sediment transport rate sequence qSu, upstream sediment transport rate time series tus, downstream hydrological station flow sequence qd, downstream flow time series tdq, downstream hydrological station sediment transport rate sequence qSd, and downstream sediment transport rate time series tud. Identify the corresponding times of the flood peaks at upstream and downstream hydrological stations during the same flood event. Determine the flood propagation time based on the time difference between the two times. The flood propagation time varies for different flood events. Select the average value to determine the flood propagation time T between the two hydrological stations. Generally, T is greater than 1 hour. Select the time interval Δt of the compiled data as an approximation of the flood propagation time or a divisor of its approximation, satisfying Δt≈T / m, where m is a positive integer. Δt is also the calculation period length in subsequent steps, generally Δt is taken as 1 hour or a multiple thereof. The upstream hydrological station flow sequence qu is expressed as:

[0097] qu={q u (t uq 1 ), q u (t uq 2 ), q u (t uq 3 ), ..., qu (t uq n1 )},

[0098] The corresponding upstream flow time series tuq is represented as:

[0099] tuq={t uq 1 , t uq 2 , t uq 3 , ..., t uq n1},

[0100] Where t uq 1 t uq 2 t uq 3 t uq n1 These represent the times 1, 2, 3, and n1 of the flow sequence, respectively. u (t uq 1 ), q u (t uq 2 ), q u (t uq 3 ), q u (t uq n1 ) represent t respectively uq 1 t uq 2 t uq 3 t uq n1 The flow rate at the upstream hydrological station at time n is n1, and the length of its flow rate sequence is n1; qu and tuq are the names of the flow rate sequence at the upstream hydrological station and the name of the upstream flow rate time series, respectively.

[0101] The sediment transport rate sequence qSu of upstream hydrological stations is represented as:

[0102] qSu={q Su (t us 1 ), q Su (t us 2 ), q Su (t us 3 ), ..., q Su (t us n2 )},

[0103] The corresponding upstream sediment transport rate time series tus is represented as:

[0104] tus={t us 1 , t us 2 , t us 3 , ..., t us n2},

[0105] Where t us 1 t us 2 t us 3 t us n2 These represent the times of the sediment transport rate sequence at the upstream hydrological stations, specifically the times of the 1st, 2nd, 3rd, and n2nd points, respectively. Su (tus 1 ), q Su (t us 2 ), q Su (t us 3 ), q Su (t us n2 ) represent t respectively us 1 t us 2 t us 3 t us n2 The sediment transport rate of upstream hydrological stations at any given time is n2, and the length of the sediment transport rate sequence of upstream hydrological stations is n2.

[0106] The downstream hydrological station flow sequence qd is represented as:

[0107] qd={q d (t dq 1 ), q d (t dq 2 ), q d (t dq 3 ), ..., q d (t dq n3 )}

[0108] The corresponding downstream flow time series tdq is represented as:

[0109] tdq={t dq 1 , t dq 2 , t dq 3 , ..., t dq n3},

[0110] Where t dq 1 t dq 2 t dq 3 t dq n3 q represents the times of the flow sequence at downstream hydrological stations at the 1st, 2nd, 3rd, and nth consecutive moments. d (t dq 1 ), q d (t dq 2 ), q d (t dq 3 ), q d (t dq n3 ) represent t respectively dq 1 t dq 2 t dq 3 t dq n3 The flow rate at the downstream hydrological station at a given time is n, and the length of the flow rate sequence at the downstream hydrological station is n3.

[0111] The sediment transport rate sequence qSd of downstream hydrological stations is represented as:

[0112] qSd={q Sd (t ds 1 ), q Sd (t ds 2), q Sd (t ds 3 ), ..., q Sd (t ds n4 )}

[0113] The corresponding downstream sediment transport rate time series tud is expressed as:

[0114] tud={t ud 1 , t ud 2 , t ud 3 , ..., t ud n4},

[0115] Where t ud 1 t ud 2 t ud 3 t ud n4 These represent the times of the sediment transport rate sequence at downstream hydrological stations, specifically the 1st, 2nd, 3rd, and nth times respectively. Sd (t ds 1 ), q Sd (t ds 2 ), q Sd (t ds 3 ), q Sd (t ds n2 ) represent t respectively ds 1 t ds 2 t ds 3 t ds n4 The sediment transport rate at downstream hydrological stations is n = 4. Here, the flow and sediment transport rate time series of upstream and downstream hydrological stations are not arithmetic sequences.

[0116] S3. Use linear interpolation to normalize the flow and sediment transport rate data sequences of upstream and downstream hydrological stations. Step S2 ensures that the data sequences from upstream and downstream hydrological stations have the same time interval, the same moment, and the same time period. Take T1 as the maximum value of the first term of the four time series for flow and sediment transport rate at upstream and downstream hydrological stations, where T1 = max{t uq1 , t us 1 , t dq 1 , t ds 1 Let T2 be the minimum value of the last term of the four time series, then T2 = min{t}. uq n1 , t us n2 , t dq n3 , t ds n4}, calculate (T2-T1) / Δt+1 and take its integer value as the length N of the unified time series; starting from t1=T1, determine the unified time series {t1, t2, t3, ..., t4} with the time interval Δt as the common difference. N}, where t1, t2, t3, t NLet represent the times 1, 2, 3, and N in sequence. This time series is an arithmetic sequence, and its start and end points are all within the range of the four time series in step 2. For any term t within the unified time series... i Let i take values ​​from 1 to N, and find t. i In step 2, at the position of the four time series, if t i In the upstream flow time series tuq, t satisfies uq k <t i <t uq k+1 Where k is some integer, then Q u (t i )=q u (t uq k )+(q u (t uq k+1 )-q u (t uq k ))(t i -t uq k ) / (t uq k+1 -t uq k ), Q u (t i ) for t i The flow rate at upstream hydrological stations at given times is used to obtain the upstream hydrological station flow rate specification sequence QU, represented as follows:

[0117] QU={Q u (t1), Q u (t2), Q u (t3), ..., Q u (t N )}

[0118] Among them, Q u (t1), Q u (t2), Q u (t3), Q u (t N ) represent t1, t2, t3, t respectively. N Similarly, the sediment transport rate standard sequence QSU for upstream hydrological stations is obtained as follows:

[0119] QSU={Q Su (t1), Q Su (t2), Q Su (t3), ..., Q Su (t N )}

[0120] Among them, Q Su (t1), Q Su (t2), Q Su (t3), Q Su (tN ) represent t1, t2, t3, t respectively. N The sediment transport rate at the upstream hydrological station is given by the following formula: Similarly, the standard sequence of flow rate QD and the standard sequence of sediment transport rate QSD at the downstream hydrological station are expressed as follows:

[0121] QD={Q d (t1), Q d (t2), Q d (t3), ..., Q d (t N )}

[0122] QSD={Q Sd (t1), Q Sd (t2), Q Sd (t3), ..., Q Sd (t N )}

[0123] Among them, Q d (t1), Q d (t2), Q d (t3), Q d (t N ) represent t1, t2, t3, t respectively. N Downstream hydrological station flow rate, Q Sd (t1), Q Sd (t2), Q Sd (t3), Q Sd (t N ) represent t1, t2, t3, t respectively. N Sediment transport rate at downstream hydrological stations at any given time;

[0124] S4. Determine the values ​​of the flow calculation coefficients using a trial-and-error method or the least squares method.

[0125] When using the trial-and-error method, the initial time t1 of the channel storage capacity WP(t1) is set to zero, and the values ​​of each term in the channel storage capacity sequence WP are calculated using the following formula:

[0126] WP={WP(t1), WP(t2), WP(t3),..., WP(t N )}

[0127] WP(t i ) = WP(t i-1 )+(Q u (t i )-Q d (t i ))Δt,i=2,3,...,N

[0128] Among them, WP(t1), WP(t2), WP(t3), WP(t N ) represent t1, t2, t3, t respectively. N The water storage capacity of the channel in the river section at any given time; WP(t) i ),WP(t) i-1 ) represent t respectively i t i-1 The water storage capacity of the river section at any given time; Q u (t i ), Q d (t i ) represent t respectively i The upstream and downstream hydrological station flows are given at specific times; Δt is the calculation period, consistent with S2. The values ​​of each term in the indicated storage flow sequence QP of the river segment are calculated using the following formula for the upstream and downstream hydrological station flow specification sequences QU and QD:

[0129] QP={QP(t1), QP(t2), QP(t3),..., QP(t N )}

[0130] QP(t i )=x Q u (t i )+(1-x)Q d (t i ), i = 1, 2, ..., N

[0131] Among them, QP(t1), QP(t2), QP(t3), QP(t) N ) represent t1, t2, t3, t respectively. N The measured flow rate of the river section at any given time; Q u (t i ), Q d (t i ) are the values ​​of the i-th item in the discharge specification sequences QU and QD of the upstream and downstream hydrological stations, respectively (i.e., t). i (Current flow at upstream hydrological stations and downstream hydrological stations), QP(t) i ) represents the value of the i-th term in the reservoir discharge sequence QP of the river segment (i.e., t). iindicates the indicated storage flow of the river reach at the time; the discharge specific gravity factor x is in the range of 0 to 1, and the set values are x1, x2, x3 respectively, satisfying 0<x1<x2<x3<1. The indicated storage flow sequences QP1, QP2, QP3 corresponding to x1, x2, x3 are calculated by using the calculation formula of indicated storage flow of the river reach. The relationship between the indicated storage flow sequences QP1, QP2, QP3 and the river reach storage capacity sequence WP is plotted respectively, and the x value that enables the closest linear relationship between the two is determined through graphic observation and comparison. When the x value is x1, x4, x5, x6 are taken, satisfying 0<x4<x5=x1<x6<x2; when the x value is x2, x4, x5, x6 are taken, satisfying x1<x4<x5=x2<x6<x3; when the x value is x3, x4, x5, x6 are taken, satisfying x2<x4<x5=x3<x6<1. The indicated storage flow sequences QP4, QP5, QP6 corresponding to x4, x5, x6 are calculated by using the calculation formula of indicated storage flow of the river reach. The relationship between the indicated storage flow sequences QP4, QP5, QP6 and the river reach storage capacity sequence WP is plotted respectively, and the x value that enables the closest linear relationship between the two is determined through graphic observation and comparison. The above process is repeated in sequence until no further comparison can be made through graphic observation. The slope value of the linear relationship between the indicated storage flow Q' and the river reach storage capacity W is the storage constant K. After K and x are determined through trial calculation, the flow routing coefficients C0, C1, C2 need to be converted and determined by the following formula:

[0132]

[0133] When calculation is performed by the least square method, the values of flow routing coefficients C0, C1, C2 can be directly solved and determined by the following matrix operation:

[0134]

[0135]

[0136] Wherein, Y is an (N-1)-dimensional column vector formed by the normalized flow sequence QD of the downstream hydrological station, and X is an (N-1)×3-dimensional matrix formed by the normalized flow sequences QD and QU of the downstream hydrological station and the upstream hydrological station, is a three-dimensional column vector formed by flow routing coefficients, X T represents the transpose matrix of X.

[0137] S5, The Muskingum method flow routing formula is used to route the flow process, and the formula is as follows:

[0138] Q pd (t i+1 )=C0Q u (t i+1 )+C1Q u (t i )+C2Q d (ti ), i = 1, 2, ..., N-1

[0139] Among them, Q u (t i+1 ), Q u (t i ) are respectively t i+1 t i Flow rate at upstream hydrological stations; Q d (t i ) for t i Downstream hydrological station flow rate; Q pd (t i+1 ) for t i+1 The calculated flow values ​​of downstream hydrological stations can form a sequence of flow calculations for downstream hydrological stations:

[0140] QPD={Q pd (t2), Q pd (t3), ..., Q pd (t N )}

[0141] The accuracy of the method is determined by comparing it with the flow specification sequence of downstream hydrological stations, according to the following formula for calculating the flow determination factor:

[0142]

[0143] in, The flow determination coefficient is the flow rate determination factor calculated using the Muskingan method flow rate calculation formula. The closer the value of Q is to 1, the higher the accuracy of the method of this invention. pd (t i ) for t i Calculated flow values ​​at downstream hydrological stations at any given time; The average value of the downstream hydrological station flow specification sequence QD, excluding the first term, is calculated using the following formula:

[0144]

[0145] S6. The sediment transport process is calculated using the Muskingen method sediment transport rate calculation formula, as follows:

[0146] Q Spd (t i+1 ) = C S,0 Q Su (t i+1 )+C S,1 Q Su (t i )+C S,2 Q Sd (t i )+CS,3 S S (t i ), i = 1, 2, ..., N-1

[0147] Among them, C S,0 ~C S,3 C is the coefficient for calculating sediment transport rate. S,0 C S,1 C S,2 These are consistent with the flow calculation coefficients C0, C1, and C2, respectively. S,3 =C S,0 +C S,1 S S (t i ) for t i The source of siltation and sedimentation in the river channel. Q Spd (t i+1 ) for t i+1 The calculated sediment transport rates of downstream hydrological stations can form a sequence of sediment transport rate calculations for downstream hydrological stations:

[0148] QSPD = {Q Spd (t2), Q Spd (t3), ..., Q Spd (t N )},

[0149] The accuracy of the method is determined by comparing it with the standard sequence of sediment transport rates at downstream hydrological stations, according to the following formula for calculating the coefficient of determination:

[0150]

[0151] in, The coefficient of determination of sediment transport rate is the formula for sediment transport rate calculation using the Muskingen method. The closer the value of Q is to 1, the higher the accuracy of the method of this invention. Spd (t i ) for t i Calculated sediment transport rate at downstream hydrological stations at any given time; The average value of the sediment transport rate standard sequence QSD for downstream hydrological stations, excluding the first term, is calculated using the following formula:

[0152]

[0153] t i Riverbed siltation source project S S (t i Based on the specific relationship between the sediment transport capacity and sediment content of the river channel, when using the unbalanced sediment transport method, we have:

[0154]

[0155] Wherein, S(t) i ), S * (t i ) represent the average suspended sediment concentration and sediment-carrying capacity of the river section, respectively; ω represents the settling velocity of sediment particles; α represents the recovery saturation coefficient of suspended sediment; h(t) represents the average suspended sediment concentration and sediment-carrying capacity of the river section, respectively; ω represents the settling velocity of sediment particles; α represents the recovery saturation coefficient of suspended sediment; h(t) represents the sediment-carrying capacity of the river section, respectively; α ... i ) represents the average water depth of the river cross section; L represents the river length between the two hydrological stations in step S1; A(t) represents the average water depth of the river cross section. i S represents the cross-sectional area of ​​the river channel. The source term can also be calculated using empirical formulas based on the actual river section conditions. When the source term is not considered, S... S (t i ) = 0.

[0156] In this method, the flow propagation time corresponds to the flood wave velocity. When calculating the sediment transport process using the Muskingan method's sediment transport rate calculation formula, the sediment is used with the same propagation time as the water flow, meaning this method assumes synchronous propagation. When their propagation velocities are inconsistent, causing the sediment peak to lag behind the flood peak, it is necessary to select the sediment propagation time for calculation. In this case, the sediment specific gravity factor x... S Still consistent with the flow proportion factor x S =x, the sediment propagation time is taken to be consistent with the time of motion of water particles, that is, the speed of sediment moving downstream is equal to the speed of water flow. The sediment propagation time is calculated using the following formula:

[0157] Ks=2L / (V um +V dm )

[0158] Among them, K S It refers to the time it takes for sediment to travel, V um V dm These are the average flow velocities at the peak flood times of the upstream and downstream hydrological stations, respectively. The following coefficients are used to calculate the sediment transport rate:

[0159]

[0160] C S,3 =C S,0 +C S,1

[0161] This invention expands upon the traditional Muskingan method flood calculation formula to develop a Muskingan sediment transport rate calculation formula. The Muskingan method sediment transport rate calculation formula for a river section is obtained by simultaneously solving the river section sediment balance equation and sediment storage equation. This invention simultaneously achieves rapid and effective calculation of water and sediment processes, has strong applicability and scalability, and can serve as a rapid and effective hydrological method for water and sediment calculation in sediment-rich rivers.

Claims

1. A river water and sediment calculation method based on the Muskingen method, characterized in that, Includes the following steps: S1. Determine the target river section to be calculated, select appropriate hydrological stations as the upstream and downstream hydrological stations of the target river section, and determine the river length L between the two hydrological stations. S2. Collect historical hydrological data from upstream and downstream hydrological stations, and compile information on flow and sediment transport rates for the same period at both stations, including the upstream hydrological station flow sequence qu, upstream flow time series tuq, upstream hydrological station sediment transport rate sequence qSu, upstream sediment transport rate time series tus, downstream hydrological station flow sequence qd, downstream flow time series tdq, downstream hydrological station sediment transport rate sequence qSd, and downstream sediment transport rate time series tud; determine the flood propagation time at the two hydrological stations, and select the time interval of the compiled data as an approximation of the flood propagation time or a divisor of its approximation; S3. Use linear interpolation to standardize the flow and sediment transport rate data sequences of upstream and downstream hydrological stations, so that the data sequences of upstream and downstream hydrological stations are at the same time interval, at the same moment, and at the same time period. S4. Use trial-and-error or least squares method to determine the values ​​of flow calculation coefficients C0, C1, and C2; S5. The flow process is calculated using the Muskingum method flow calculation formula to obtain the flow calculation sequence of downstream hydrological stations; S6. The sediment transport process was calculated using the Muskingan method sediment transport rate calculation formula to obtain the sediment transport rate calculation sequence for downstream hydrological stations; the Muskingan method sediment transport rate calculation formula is as follows: Q Spd (t i+1 )=C S,0 Q Su (t i+1 )+C S,1 Q Su (t i )+C S,2 Q Sd (t i )+C S,3 S S (t i ),i=1,2,…,N-1 Among them, C S,0 ~C S,3 C is the coefficient for calculating sediment transport rate. S,0 C S,1 C S,2 The values ​​of C0, C1, and C2 are the same as those of the flow calculation coefficients, respectively. S,3 =C S,0 +C S,1 S S (t i ) for t i Riverbed erosion and sedimentation source; Q Su (t i+1 ), Q Su (t i ) are respectively t i+1 t i Sediment transport rate at upstream hydrological stations at any given time; Q Sd (t i ) for t i Sediment transport rate at downstream hydrological stations at any given time; Q Spd (t i+1 ) for t i+1 The calculated sediment transport rate values ​​at downstream hydrological stations at given time points constitute the sediment transport rate calculation sequence for downstream hydrological stations as follows: QSPD={Q Spd (t2), Q Spd (t3), ..., Q Spd (t N )}; The accuracy of the method is determined by comparing it with the standard sequence of sediment transport rates at downstream hydrological stations, according to the following formula for calculating the coefficient of determination: in, The coefficient of determination of sediment transport rate is the formula for sediment transport rate calculation using the Muskingen method. The closer the value of Q is to 1, the higher the calculation accuracy. Spd (t i ) for t i Calculated sediment transport rate at downstream hydrological stations at any given time; This is the average value of the sediment transport rate standard sequence QSD for downstream hydrological stations, excluding the first term.

2. The river water and sediment calculation method based on the Muskingan method according to claim 1, characterized in that, Step S2, which involves compiling the flow and sediment transport rate information of the two hydrological stations during the same period, specifically includes: The upstream hydrological station flow sequence qu is represented as: qu={q u (t uq 1 ), q u (t uq 2 ), q u (t uq 3 ), ..., q u (t uq n1 )}; The corresponding upstream flow time series tuq is represented as: tuq={t uq 1 , t uq 2 , t uq 3 , ..., t uq n1 }; Among them, t uq 1 t uq 2 t uq 3 t uq n1 These represent the times of the flow sequence at the upstream hydrological stations at the 1st, 2nd, 3rd, and nth consecutive moments, respectively. u (t uq 1 ), q u (t uq 2 ), q u (t uq 3 ), q u (t uq n1 ) represent t respectively uq 1 t uq 2 t uq 3 t uq n1 The flow rate at the upstream hydrological station at a given time, with a flow sequence length of n1; The sediment transport rate sequence qSu of upstream hydrological stations is represented as: qSu={q Su (t us 1 ), q Su (t us 2 ), q Su (t us 3 ), ..., q Su (t us n2 )}; The corresponding upstream sediment transport rate time series tus is represented as: tus={t us 1 , t us 2 , t us 3 , ..., t us n2 }; Among them, t us 1 t us 2 t us 3 t us n2 These represent the times of the sediment transport rate sequence at the upstream hydrological stations, specifically the times of the 1st, 2nd, 3rd, and n2nd points, respectively. Su (t us 1 ), q Su (t us 2 ), q Su (t us 3 ), q Su (t us n2 ) represent t respectively us 1 t us 2 t us 3 t us n2 The sediment transport rate of upstream hydrological stations at any given time, and the length of the sediment transport rate sequence of upstream hydrological stations is n². The downstream hydrological station flow sequence qd is expressed as: qd={q d (t dq 1 ), q d (t dq 2 ), q d (t dq 3 ), ..., q d (t dq n3 )}; The corresponding downstream flow time series tdq is represented as: tdq={t dq 1 , t dq 2 , t dq 3 , ..., t dq n3 }; Among them, t dq 1 t dq 2 t dq 3 t dq n3 q represents the times of the flow sequence at downstream hydrological stations at the 1st, 2nd, 3rd, and nth consecutive moments. d (t dq 1 ), q d (t dq 2 ), q d (t dq 3 ), q d (t dq n3 ) represent t respectively dq 1 t dq 2 t dq 3 t dq n3 The flow rate at downstream hydrological stations at given time, and the length of the flow rate sequence at downstream hydrological stations is n3; The sediment transport rate sequence qSd at downstream hydrological stations is expressed as: qSd={q Sd (t ds 1 ), q Sd (t ds 2 ), q Sd (t ds 3 ), ..., q Sd (t ds n4 )}; The corresponding downstream sediment transport rate time series tud is expressed as: tud={t ud 1 , t ud 2 , t ud 3 , ..., t ud n4 }; Among them, t ud 1 t ud 2 t ud 3 t ud n4 These represent the times of the sediment transport rate sequence at downstream hydrological stations, specifically the 1st, 2nd, 3rd, and nth times respectively. Sd (t ds 1 ), q Sd (t ds 2 ), q Sd (t ds 3 ), q Sd (t ds n2 ) represent t respectively ds 1 t ds 2 t ds 3 t ds n4 The sediment transport rate of downstream hydrological stations at a given time, and the length of the sediment transport rate sequence of downstream hydrological stations is n4; The method for calculating the flood propagation time T is as follows: Based on the compiled information on flow and sediment transport rate of the two hydrological stations during the same period, the corresponding times of the flood peaks at the upstream and downstream hydrological stations during the same flood event are identified. The time difference between the two is used to determine the flood propagation time of the event. The flood propagation time is different for different flood events. The average value is selected to determine the flood propagation time T between the two hydrological stations. The time interval Δt of the compiled data is chosen as an approximation of the flood propagation time or a divisor of that approximation.

3. The river water and sediment calculation method based on the Muskingan method according to claim 1, characterized in that, Step S3 is as follows: Let T1 be the maximum value of the first term of the four time series of flow and sediment transport rate at upstream and downstream hydrological stations, T1 = max{t uq 1 , t us 1 , t dq 1 , t ds 1 Let T2 be the minimum value of the last term of the four time series, then T2 = min{t}. uq n1 , t us n2 , t dq n3 , t ds n4 }, calculate (T2-T1) / Δt+1 and take its integer value as the length N of the unified time series; starting from t1=T1, determine the unified time series {t1, t2, t3, ..., t4} with the time interval Δt as the common difference. N }, where t1, t2, t3, t N Let t represent the times 1, 2, 3, and N in sequence. This time series is an arithmetic sequence, and its start and end points are all within the four time series ranges of flow and sediment transport rate at upstream and downstream hydrological stations. For any term t within the unified time series... i Let i take values ​​from 1 to N, and find t. i The positions of the flow and sediment transport rate at the upstream and downstream hydrological stations in four time series, if t i In the upstream flow time series tuq, t satisfies uq k <t i <t uq k+1 Where k is some integer, then Q u (t i )=q u (t uq k )+(q u (t uq k+1 )-q u (t uq k ))(t i -t uq k ) / (t uq k+1 -t uq k ), Q u (t i ) for t i The flow rate at upstream hydrological stations is used to obtain the upstream hydrological station flow rate specification sequence QU, which is represented as: QU={Q u (t1), Q u (t2), Q u (t3), ..., Q u (t N )} Among them, Q u (t1), Q u (t2), Q u (t3), Q u (t N ) represent t1, t2, t3, t respectively. N Similarly, the sediment transport rate standard sequence QSU for upstream hydrological stations is obtained as follows: QSU={Q Su (t1),Q Su (t2),Q Su (t3),……,Q Su (t N )} Among them, Q Su (t1), Q Su (t2), Q Su (t3), Q Su (t N ) represent t1, t2, t3, t respectively. N Sediment transport rate at upstream hydrological stations at all times; Similarly, the downstream hydrological station flow specification sequence QD and sediment transport rate specification sequence QSD are respectively expressed as: QD={Q d (t1),Q d (t2),Q d (t3),……,Q d (t N )} QSD={Q Sd (t1),Q Sd (t2),Q Sd (t3),……,Q Sd (t N )} Among them, Q d (t1), Q d (t2), Q d (t3), Q d (t N ) represent t1, t2, t3, t respectively. N Downstream hydrological station flow rate, Q Sd (t1), Q Sd (t2), Q Sd (t3), Q Sd (t N ) represent t1, t2, t3, t respectively. N Sediment transport rate at downstream hydrological stations at any given time.

4. The river water and sediment calculation method based on the Muskingan method according to claim 1, characterized in that, In step S4, when using the trial-and-error method, the initial water storage capacity WP(t1) of the river segment at time t1 is set to zero. The values ​​of each term in the water storage capacity sequence WP of the river segment are calculated using the following formula: WP={WP(t1), WP(t2), WP(t3), ..., WP(t... N )},WP(t i ) = WP(t i-1 )+(Q u (t i )-Q d (t i ))Δt, i=2, 3,..., N; Among them, WP(t1), WP(t2), WP(t3), WP(t N ) represent t1, t2, t3, t respectively. N The water storage capacity of the channel in the river section at any given time; WP(t) i ),WP(t) i-1 ) represent t respectively i t i-1 The water storage capacity of the river section at any given time; Q u (t i ), Q d (t i ) represent t respectively i The flow rates at upstream and downstream hydrological stations at specific times; the time interval Δt is the length of the calculation period. The values ​​of each term in the reservoir discharge sequence QP of the river segment are calculated using the following formula for the discharge specification sequences QU and QD of the upstream and downstream hydrological stations: QP={QP(t1), QP(t2), QP(t3), …, QP(t… N )},QP(t i )=x Q u (t i )+(1-x)Q d (t i ), i = 1, 2, ..., N; Among them, QP(t1), QP(t2), QP(t3), QP(t) N ) represent t1, t2, t3, t respectively. N The storage flow rate of the river section at any given time; Q u (t i ), Q d (t i ) are respectively t i The flow rates at upstream and downstream hydrological stations at specific times, QP(t) i ) for t i The storage flow rate of the river section at any given time; The flow specific gravity factor x is within the range of 0 to 1, and the set values are x₁, x₂, x₃ respectively, satisfying 0 < x₁ < x₂ < x₃ < 1. The indicator flow sequences QP₁, QP₂, QP₃ corresponding to x₁, x₂, x₃ are calculated by using the indicator flow calculation formula of the river reach. The relationships between the indicator flow sequences QP₁, QP₂, QP₃ and the river reach channel storage volume sequence WP are plotted respectively, and the x value that enables the two to have the closest linear relationship is determined through observation and comparison of the graphs. When the x value is x₁, x₄, x₅, x₆ are taken, satisfying 0 < x₄ < x₅ = x₁ < x₆ < x₂; when the x value is x₂, x₄, x₅, x₆ are taken, satisfying x₁ < x₄ < x₅ = x₂ < x₆ < x₃; when the x value is x₃, x₄, x₅, x₆ are taken, satisfying x₂ < x₄ < x₅ = x₃ < x₆ < 1. The indicator flow sequences QP₄, QP₅, QP₆ corresponding to x₄, x₅, x₆ are calculated by using the indicator flow calculation formula of the river reach. The relationships between the indicator flow sequences QP₄, QP₅, QP₆ and the river reach channel storage volume sequence WP are plotted respectively, and the x value that enables the two to have the closest linear relationship is determined through observation and comparison of the graphs. The above process is repeated in sequence until no further comparison can be made by observing the graphs; the slope value of the linear relationship between the indicator flow Q' and the river reach channel storage volume W is the storage constant K. After K and x are determined by trial calculation, the flow routing coefficients C₀, C₁, C₂ need to be converted and determined by the following formula:

5. The river water and sediment calculation method based on the Muskingen method according to claim 1, characterized in that, When the least square method is used for calculation in step S4, the following matrix operation is directly used to solve and determine the values of flow routing coefficients C₀, C₁, C₂: Where Y is an N-1 dimensional column vector composed of the downstream hydrological station flow specification sequences QD, and X is an (N-1)×3 dimensional matrix composed of the downstream and upstream hydrological station flow specification sequences QD and QU. X is a three-dimensional column vector composed of flow calculation coefficients. T Let X be the rank matrix of X.

6. The river water and sediment calculation method based on the Muskingan method according to claim 1, characterized in that, In step S5, the Muskingum flow routing formula is: Q pd (t i+1 )=C0Q u (t i+1 )+C1Q u (t i )+C2Q d (t i ),i=1,2,…,N-1 Among them, Q u (t i+1 ), Q u (t i ) are respectively t i+1 t i Flow rate at upstream hydrological stations at all times; Q d (t i ) for t i Downstream hydrological station flow rate; Q pd (t i+1 ) for t i+1 The calculated flow values ​​at downstream hydrological stations can form a flow calculation sequence for downstream hydrological stations: QPD={Q pd (t2), Q pd (t3), ..., Q pd (t N )}; The accuracy of the method is determined by comparing it with the flow specification sequence of downstream hydrological stations, according to the following formula for calculating the flow determination factor: in, Q is the flow determination factor calculated using the Muskingan method flow calculus formula. pd (t i ) for t i Calculated flow values ​​at downstream hydrological stations at any given time; The flow specification sequence QD for downstream hydrological stations is the average value excluding the first term; the flow determination coefficient. The closer the value is to 1, the higher the calculation accuracy.

7. The river water and sediment calculation method based on the Muskingan method according to claim 1, characterized in that, In step S6, the source term is determined according to the relationship between the sediment transport capacity and sediment concentration of the specific river channel. When the unbalanced sediment transport form is adopted, there is: Wherein, S(t) i ), S * (t i ) represent the average suspended sediment concentration and sediment-carrying capacity of the river section, respectively; ω represents the settling velocity of sediment particles; α represents the recovery saturation coefficient of suspended sediment; h(t) represents ... i ) represents the average water depth of the river cross section; L represents the river length between the two hydrological stations; A(t) represents the average water depth of the river cross section; L represents the average water depth of the river cross section between the two hydrological stations; A(t) represents the average water depth of the river cross section between the two hydrological stations. i () represents the cross-sectional area of ​​the river channel. Alternatively, the source term can be calculated using empirical formulas based on the actual river section conditions. When the source term is not considered, S S (t i =0, and the formula for calculating sediment transport rate and flow rate is consistent.

8. A river water and sediment calculation system based on the Muskingen method, characterized in that, comprising: a data acquisition and processing module, configured to determine a target river reach, select an upstream hydrological station and a downstream hydrological station of the target river reach, and determine the river length between the two hydrological stations; collecting the historical hydrological data of the upstream hydrological station and the downstream hydrological station, compiling the flow and sediment concentration information of the two hydrological stations in the same period, determining the flood propagation time between the two stations, and selecting the time interval of the compiled data as an approximate value or a divisor of the flood propagation time; compiling the flow and sediment concentration (or sediment discharge) data sequences of the upstream hydrological station and the downstream hydrological station, and performing linear interpolation on the measured data so that the data sequences satisfy that the data of the upstream and downstream hydrological stations have the same time interval, the same time point and the same time period; a flow routing coefficient calculation module, configured to determine the values of flow routing coefficients by a trial algorithm or the least square method; a flow process routing module, configured to route the flow process by using the Muskingum flow routing formula to obtain the calculated value of the downstream station flow process; a sediment transport process routing module, configured to route the sediment transport process by using the Muskingum sediment discharge routing formula to obtain the calculated value of the downstream station sediment transport process.

9. A device, characterized in that, comprising a memory and a processor, wherein: the memory is configured to store a computer program capable of running on the processor; the processor is configured to, when running the computer program, execute the steps of the river water and sediment routing method based on the Muskingum method according to any one of claims 1 to 7.

10. A storage medium, characterized in that, The storage medium stores a computer program that, when executed by at least one processor, implements the steps of the river water and sediment calculation method based on the Muskingan method as described in any one of claims 1-7.

Images