Listen

In this paper, the augmented version of finite volume Harten–Lax–van Leer (HLL) solver, including source terms, is extended to free-surface and pressurized mixed pipe flows over complex and frictional topography. This augmented HLL Riemann solver is employed for the flux approximation at the cell interface, where source terms are split into two parts based on the wave propagation speed. The friction term is treated using a splitting implicit method to maintain stability over dry beds. The Preissmann slot method is adopted to reproduce pressurized flow in the conduit. The performances of the numerical model are investigated by several numerical tests and compared with existing methods showing clear improvements.

Listen

  • An augmented HLL Riemann solver is extended to mixed pipe flows over complex topography.

  • Superiority simulating steady flows over complex topography with different hydraulic flow regimes.

  • Capability simulating transient dambreak pipe flow, including wetting or drying front and the calculation of the friction term.

augmented HLL Riemann solver, mixed pipe flow, Preissmann slot method, shallow water equations
Listen

In the context of climate change and rapid urbanization over the last decades, risks of urban flooding have greatly increased worldwide with increased runoff at local catchment, that devastating consequences occurred on both infrastructures and human safety (Zheng et al. 2015). To improve risk assessment and evacuation management, flood inundation models are crucial tools to predict flood hydrodynamics in cities (Mignot et al. 2006; Hunter et al. 2007; Arrault et al. 2016; Teng et al. 2017; Chen et al. 2018). As an integrated component of urban flooding processes, drainage networks composed of junctions and pipes are commonly modeled by different numerical approaches (Martins et al. 2017; Li et al. 2020). Pipe flows, which are the basic element of drainage network, are often modeled by the 1D shallow water equations (Kerger et al. 2011a, 2011b).

During urban flooding, transitions in pipes between free-surface and pressurized flow conditions can take place repeatedly. This phenomenon has been studied by many researchers (Bousso et al. 2012) and can be modeled using a method like the shock-fitting method (Wiggert 1972) or the shock-capturing method (Abbott & Cunge 1982; Vasconcelos et al. 2006; Vasconcelos et al. 2009). The extension of free-surface shallow water equations to the pressurized flow to capture the interface has made the shock-capturing method more attractive than the cumbersome shock-fitting method. Among them, the Preissmann slot method (PSM) has been widely used for its simplicity in model implementation. Lately, Kerger et al. (2011a, 2011b) extended the classical PSM by developing a negative slot to reproduce subatmospheric pressurized flow, and Maranzoni et al. (2015) extended 1D PSM to 2D accounting for 2D pressurized flow under the bridge. Except the PSM, Vasconcelos et al. (2006) and Vasconcelos et al. (2009) developed a two-component pressure approach with the assumption that the pipe walls are elastic. In this way, various unsteady flows can be simulated effectively, including free-surface flow, pressurized flow, mixed flow and sub-atmospheric pressurized flow. PSM is used in the current study for mixed pressurized and free-surface pipe flows.

1D shallow water equations may be solved numerically by using a Godunov-type finite volume scheme (Toro & Garcia-Navarro 2007). When tackling complex topography in the framework of Godunov-type scheme, the slope source term needs to be treated properly to reproduce flow hydrodynamics. A number of numerical techniques dealing with the source term have been reviewed in the literature (Guinot 2012) to construct well-balanced schemes. For instance, Liang & Marche (2009) proposed a pre-balanced shallow water equation that no special treatment is needed for the source term considering pressure balancing. For circular cross-sections like pipe flow, Capart et al. (2003) reconstructed the momentum flux to account for slope and non-prismaticity by considering the balance of hydrostatic pressure with the approximated water surface level to achieve hydrostatic equilibrium, which was also employed for mixed pipe flow modeling by Sanders & Bradford (2010) and Aureli et al. (2015) over complex topography. The balanced property of the above-mentioned schemes is achieved by forcing equilibrium among fluxes and source terms in cases of quiescent water. However, errors on computed discharge will accumulate around the bottom variation, which can be observed in the results published in Liang & Borthwick (2009). To overcome this problem, Murillo & García-Navarro (2012) presented augmented versions of the HLL Riemann solvers for shallow water flows, by incorporating an extra wave that accounts for the presence of the source term and improved results have been obtained. In the work of Murillo & García-Navarro (2010), LeVeque (2011) and Murillo & García-Navarro (2012), only rectangular cross sections in 1D open channel flows are considered.

Those previous works are extended in the present contribution to circular cross-sections and mixed pipe flows over complex topography. Based on Murillo & García-Navarro (2010) and Murillo & García-Navarro (2012), a source term adapted to circular sections has been added to the continuity equation that along with the numerical scheme, and upgraded celerity calculations, leads to improved simulations of free-surface and mixed circular pipe flows. Following Murillo & García-Navarro (2010) and Murillo & García-Navarro (2012), the source term is discretized approximately as a succession of discontinuities separating regions of constant state and split into two parts based on the wave propagation speed. Governing equations are solved using a Godunov-type scheme with an augmented HLL approximate Riemann solve developed for mixed flow (León et al. 2006). The friction term is included and discretized by a splitting implicit method (Liang & Marche 2009) to maintain stability for flows over dry beds.

The paper is organized as follows: 1D shallow water equations with conserved variables A (flow area) and Q (total discharge) are presented in the ‘Governing equations’ section with the introduction of PSM for mixed flow modeling. In the ‘Numerical scheme’ section, a first-order Godunov-type scheme with augmented HLL Riemann solver, including source terms, is extended to the studied governing system for mixed pipe flows. The adaptations of the numerical method presented for simulating free-surface and mixed circular pipe flows are tested against the method from Sanders & Bradford (2010) in the ‘Numerical experiments’ section. Conclusions are drawn in the ‘Conclusions’ section.

Listen
A hyperbolic conservation law of the 1D shallow water equations with the conserved variables A and Q may be expressed in a matrix form as:
(1)
where t denotes time and x is the distance. ⁠, and are the vectors of conserved variables, fluxes and source terms, respectively, and are given by:
(2)
where is the flow area, is the total discharge, is the gravitational acceleration, is the bottom slope term (⁠ is the bottom elevation) and I is the hydrostatic pressure force term, which is defined as:
(3)
where h represents the water depth, and is the cross-sectional width at elevation above the bottom (Guinot 2012). is the friction slope and may evaluated as:
(4)
where is Manning's roughness coefficient, and is the hydraulic radius.
To simulate mixed flow with the classical shallow water equations in closed conduit, the PSM is used that a narrow slot is inserted into the crown of the conduit. The slot width can be determined based on the acoustic wave velocity, which is defined as:
(5)
where is the slot width, is the conduit cross-sectional area, and a is the conduit acoustic wave velocity. For circular pipe with diameter ⁠, the flow area A and pressure term I for free-surface and pressurized flow based on PSM are given, respectively, as follows (Sturm 2010):
(6)
(7)
where is the wetted angle.
Listen
In a finite volume first-order Godunov-type scheme, the variables can be updated to the next time step based on the following formula:
(8)
where the subscript i is the cell index, and the superscript n is the time index. and are the time step and cell size, respectively. and are the fluxes calculated at the interface and ⁠, respectively, and is the source term.

Numerical flux calculation

Listen
Numerical fluxes are estimated by solving the Riemann problem at the cell interface. For the 1D mixed flow modeling, interface fluxes and are calculated with the HLL approximate Riemann solver (Harten et al. 1983; León et al. 2006; Sanders & Bradford 2010). Taking the cell interface as an example, the fluxes can be evaluated by solving a local Riemann problem, i.e., ⁠, in the context of a Godunov-type scheme. and are the left and right Riemann states, which is used to define the Riemann problem and hence calculate the fluxes. The solution of the Riemann problem is approximated by an intermediate region of constant state separated from the left and right states and by two shocks. Therefore, the numerical flux is (Toro 2001):
(9)
where and are the left and right characteristic wave speeds calculated by (Toro 1992; Fraccarollo & Toro 1995):
(10)
in which u is the averaged flow velocity defined as and is calculated by:
(11)
For the circular pipe, is related to and its approximations are given by Sanders & Bradford (2010):
(12)
and subsequently,
(13)
When ⁠, the flow is under free-surface condition and the intermediate wave speed is calculated using
(14)
with
(15)
When ⁠, the water surface level may reach the crown of the pipe and the flow becomes pressurized. The intermediate wave speed is then set to be the speed of the gravity wave c calculated based on the slot width ⁠.

Bottom slope source term calculation

Listen
The stationary-preserving method to treat the bottom slope source term used here is first presented by Capart et al. (2003) for natural channel and further developed by Sanders & Bradford (2010) for mixed flow modeling. As illustrated in Figure 1, by assuming the slope of the water surface is mild and the streamwise component of the sidewall inclination is small, momentum balance is approximated using an average elevation of water surface within the control volume instead of the local free-surface elevation ⁠. The overall pressure thrust acting on control volume is approximated as:
(16)
Figure 1
Hydrostatic balance of channel reach featuring irregular bathymetry adapted from Capart et al. (2003), where is the free-surface elevation.
View largeDownload slide

Hydrostatic balance of channel reach featuring irregular bathymetry adapted from Capart et al. (2003), where is the free-surface elevation.

Figure 1
Hydrostatic balance of channel reach featuring irregular bathymetry adapted from Capart et al. (2003), where is the free-surface elevation.
View largeDownload slide

Hydrostatic balance of channel reach featuring irregular bathymetry adapted from Capart et al. (2003), where is the free-surface elevation.

Close modal
As presented by Capart et al. (2003), when the channel is rectangular in shape and has a constant width (regardless of the slope of the water surface) and when the water surface is horizontal (regardless of channel geometry), momentum balance is achieved exactly. The bottom slope source term is treated explicitly and is calculated as:
(17)
where ⁠.

In the following, the above-described schemes with Equation (17) are referred to as the ‘Sanders & Bradford (2010)’ scheme.

Augmented HLL Riemann solver including slope source term

Listen

In this section, the numerical scheme on source term discretization presented by Murillo & García-Navarro (2010) and Murillo & García-Navarro (2012) is extended from rectangular shape to the currently studied pipe flow system, which is based on the projection of the source term onto the Jacobian's eigenvector basis.

The introduction of source term in the Riemann solver has added a stationary shock wave into the source at the discontinuity point Two inner constant states and are introduced as in Figure 2. Therefore, the selection of estimation for the wave speed cannot be directly applied to the cases with the source term. As suggested by Murillo & García-Navarro (2012), the wave speeds and are approximated using:
(18)
(19)
where ⁠, and are wave estimates of Roe (1981) in Riemann problems when there are source terms (variable bottom topography), which are calculated as:
(20)
Figure 2
Augmented HLL Riemann solutions consisted of two inner constant states separated by a stationary shock wave at .
View largeDownload slide

Augmented HLL Riemann solutions consisted of two inner constant states separated by a stationary shock wave at ⁠.

Figure 2
Augmented HLL Riemann solutions consisted of two inner constant states separated by a stationary shock wave at .
View largeDownload slide

Augmented HLL Riemann solutions consisted of two inner constant states separated by a stationary shock wave at ⁠.

Close modal
As the approximate Riemann solution consists of only discontinuities and is constructed as a sum of jumps or shock for Roe's linearization, the Harten-Hyman entropy fix (Toro 2013) is used to avoid unphysical results. In the case of left transonic rarefaction ⁠, where and ⁠, the wave speed can be approximated as:
(21)
For a right transonic rarefraction ⁠, where and ⁠, the wave speed can be approximated as:
(22)
The slope source term is discretized approximately as a succession of discontinuities separating regions of constant state as the space derivative of the flux in upwinding schemes (Guinot 2012):
(23)
where is a vector formed by components of the wave strengths of the source term, and is the eigenvectors of 1D shallow water Equation (1) with the characteristic form with ⁠:
(24)

The total source term in Equation (23) is split into two parts:

  • one part that corresponding to negative wave speeds, ⁠, is assigned to the cell on the left-hand side of the initial discontinuity.

  • another part that corresponding to positive wave speeds, ⁠, is assigned to the cell on the right-hand side of the initial discontinuity.

Therefore, the source term can be split into two parts based on the wave speeds:
(25)
The discretization of Equation (8) can be rewritten as:
(26)
The derivation starts from the calculation of wave strengths at interface between the cell and ⁠:
(27)
The source term at the interface is split into two terms and based on the wave speeds:
(28)
(29)
with
(30)
(31)

In subcritical flow, the wave propagates to both the upstream and the downstream. Therefore, and ⁠.

For the momentum term, the second component of the flux at the interface is,
(32)
The flow area at the interface between the cell and i is reset to ⁠:
(33)
The same derivation can be applied to ⁠, and we can get,
(34)
For the continuity equation, the source terms developed by Murillo & García-Navarro (2012) are added to the Riemann solver as:
(35)
(36)

Three flow conditions related to dry bed are illustrated in Figure 3:

  • ⁠, and

Figure 3
Three flow conditions with dry bed at cell interface . is the minimum water depth for dry bed.
View largeDownload slide

Three flow conditions with dry bed at cell interface ⁠. is the minimum water depth for dry bed.

Figure 3
Three flow conditions with dry bed at cell interface . is the minimum water depth for dry bed.
View largeDownload slide

Three flow conditions with dry bed at cell interface ⁠. is the minimum water depth for dry bed.

Close modal
As illustrated in Figure 3(a), the left cell is dry and the right cell is wet but not high enough to flood the left cell. In this flow condition, no flux exchange takes place. To balance the flux gradient and source term, the source terms are given as:
(37)

  • ⁠, and

As illustrated in Figure 3(b), for the flow condition that right cell is dry and left cell is wet and no flux exchange occurs, the source terms are given as:
(38)

  • and

As illustrated in Figure 3(c), for the flow conditions that both cells are dry, all the terms are calculated as as below:
(39)

Friction calculation

Listen
The implicit treatment for friction term is necessary to guarantee the stable simulation in the very shallow cells due to that such terms become stiff with very small water depth in its denominator. The friction source terms are evaluated using a splitting point-implicit scheme for better stability followed by Liang & Marche (2009), which is equivalent to solve the following ordinary differential equations:
(40)
where with ⁠. The friction term can be isolated and expanded for the momentum equation using a Taylor series as:
(41)
where ⁠. Ignoring the higher-order terms, the following time-marching formula can then be obtained for ⁠:
(42)
where implicit coefficient, and is the friction source term including the implicit coefficient. By noting that the maximum effect of friction force is to stop the fluid, i.e., the criteria for can be derived (Liang & Marche 2009):
(43)

The value of will be replaced by the critical one when it is computed beyond the limit.

Listen

In this section, the performance of the proposed numerical model is tested on a number of numerical examples with a circular cross-section. Pipe flow is simulated in steady or transient state over complex topography. The numerical results obtained by using the ‘Sanders & Bradford (2010)’ scheme on a very fine grid are adopted as references for model comparisons. The stability of the numerical scheme is governed by the Courant–Friedrichs–Lewy criterion as the scheme is explicit. For all the test cases considered in this work, fixed timestep is used and the Courant number is less than 0.9 with and ⁠.

‘Lake at rest’ steady states

Listen
In this numerical test, the ability of the proposed numerical model to preserve the ‘Lake at rest’ steady state for irregular cross-sections is verified, like pipe flow. To this end, a pipe flow over irregular bottom topography with one bump is thus designed, which is defined as:
(44)

The initial conditions are defined as ⁠, over the whole computational domain so that the bump is completely submerged. The computed hydraulic variables (⁠ and ⁠) at are shown in Figure 4. Results show that both methods can preserve ‘lake at rest’ steady states for water surface elevation, but the proposed model outperforms the ‘Sanders & Bradford (2010)’ model considering the total discharge computed.

Figure 4
Computed results at ‘Lake at rest’ steady state.
View largeDownload slide

Computed results at ‘Lake at rest’ steady state.

Figure 4
Computed results at ‘Lake at rest’ steady state.
View largeDownload slide

Computed results at ‘Lake at rest’ steady state.

Close modal

Free-surface pipe flow over a bump

Listen
This test case aims to evaluate the performance of the proposed model to deal with a bump for free-surface pipe flow in steady state. A permanent flow is thus simulated with a spatially varying bottom elevation and without friction. The ability of the numerical scheme to deal with variable bottom is tested for a subcritical (T2a), transcritical (T2b) and hydraulic jump (T2c). The bed elevation is defined as Equation (45) with pipe diameter ⁠:
(45)

The reference solution for the test case is obtained on a very fine grid using the ‘Sanders & Bradford (2010)’ model. The reference for test case T2a (subcritical flow) is established using and ⁠. For the hydraulic jump (T2c), the discharge and ⁠. For the transcritical case (T2b), the discharge ⁠, and no water depth is necessary at the boundary since the flow is controlled by a critical depth (such as ⁠) that appears at (middle of the bump).

For three hydraulic regimes, the domain length is with a mesh of uniform cell ⁠. The initial velocity of this test is set to be zero. An initial water surface elevation ζini = 0.4 m is set for test case T2a and T2c, and ζini = 0.2 m is given for test case T2b. The computational parameters are given in Table 1. The simulation is carried out until steady flow conditions are reached over the computational domain (after barely ⁠).

Table 1

Parameters for the test case over a bump (T2)

SymbolMeaningValue
 Pipe diameter  
 Domain length  
 Bump length  
 Middle of bump  
 Minimum bottom elevation  
 Maximum bottom elevation  
 Cell size  
 Initial water surface elevation (T2a)  
 Initial water surface elevation (T2b)  
 Initial water surface elevation (T2c)  
 Prescribed total discharge at the upstream boundary   
 For subcritical flow (T2a)  
 For transcritical flow (T2b)  
 For transcritical flow with a jump (T2c)  
 Prescribed water depth at the downstream boundary   
 For subcritical flow (T2a)  
 For transcritical flow (T2b) Free-outlet 
 For transcritical flow with a jump (T2c)  
SymbolMeaningValue
 Pipe diameter  
 Domain length  
 Bump length  
 Middle of bump  
 Minimum bottom elevation  
 Maximum bottom elevation  
 Cell size  
 Initial water surface elevation (T2a)  
 Initial water surface elevation (T2b)  
 Initial water surface elevation (T2c)  
 Prescribed total discharge at the upstream boundary   
 For subcritical flow (T2a)  
 For transcritical flow (T2b)  
 For transcritical flow with a jump (T2c)  
 Prescribed water depth at the downstream boundary   
 For subcritical flow (T2a)  
 For transcritical flow (T2b) Free-outlet 
 For transcritical flow with a jump (T2c)  
View Large

The simulation results of three flow regimes are presented in Figure 5. Comparing with the reference, numerical results show that the proposed model provides more accurate predictions of the computed hydraulic variables (⁠ and ⁠) than the ‘Sanders & Bradford (2010)’ model, whatever the hydraulic flow regime is. For subcritical flow, the proposed method provides better estimations of the water surface elevation above and upstream of the bump, which can be clearly observed in the zoomed comparison in rectangular. For transcritical and hydraulic jump, the proposed method has better estimation of water surface elevation upstream the bump and more accurate shock position for hydraulic jump. Regarding the water discharge, the proposed method provides the exact value for sub- and transcritical flow, except one cell for hydraulic jump, which is located downstream the top of the bump. For the ‘Sanders & Bradford (2010)’ model, the errors accumulated upstream and downstream of the top point are due to the HLL solver adopted, which can also be observed in the results presented by Liang & Borthwick (2009).

Figure 5
Comparison of simulated water surface elevation (left), zoomed comparison in rectangular (middle) and total discharge (right) of free-surface flow over a bump: subcritical flow (top), transcritical flow (middle) and hydraulic jump (bottom).
View largeDownload slide

Comparison of simulated water surface elevation (left), zoomed comparison in rectangular (middle) and total discharge (right) of free-surface flow over a bump: subcritical flow (top), transcritical flow (middle) and hydraulic jump (bottom).

Figure 5
Comparison of simulated water surface elevation (left), zoomed comparison in rectangular (middle) and total discharge (right) of free-surface flow over a bump: subcritical flow (top), transcritical flow (middle) and hydraulic jump (bottom).
View largeDownload slide

Comparison of simulated water surface elevation (left), zoomed comparison in rectangular (middle) and total discharge (right) of free-surface flow over a bump: subcritical flow (top), transcritical flow (middle) and hydraulic jump (bottom).

Close modal

Mixed pipe flow over a bump

Listen
In this test case, mixed pipe flow modeling is performed to evaluate the performance of the proposed model to deal with a bump in steady state, and a permanent flow is thus simulated with a spatially varying bottom elevation without friction. The maximum bottom elevation has been increased to for the mixed flow, and the system is assumed to be ventilated. Two flow regimes are tested for free-surface condition, which are subcritical flow (T3a) and hydraulic jump (T3b), and one for pressurized flow (T3c). The bed elevation is defined as Equation (46) with pipe diameter and ⁠:
(46)

The reference solution for this test case is obtained on a very fine grid using the ‘Sanders & Bradford (2010)’ model. The reference for test case T3a (subcritical flow) is established using and ⁠. For the hydraulic jump (T3b), the discharge and ⁠. For the pressurized flow, ⁠, and free-outlet is used at the right side of the boundary. Same numerical parameters as free-surface pipe flow are used for this test case, which are listed in Table 2. The simulation is carried out until steady flow conditions are reached over the computational domain (after barely ⁠).

Table 2

Parameters for the test case over a bump (T3)

SymbolMeaningValue
 Pipe diameter  
 Preissmann slot width  
 Domain length  
 Bump length  
 Middle of bump  
 Minimum bottom elevation  
 Maximum bottom elevation  
 Cell size  
 Initial water surface elevation  
 Prescribed total discharge at the upstream boundary   
 For subcritical flow (T3a)  
 For hydraulic jump (T3b)  
 For pressurized flow (T3c)  
 Prescribed water depth at the downstream boundary   
 For subcritical flow (T3a)  
 For hydraulic jump (T3b)  
 For pressurized flow (T3c) Free-outlet 
SymbolMeaningValue
 Pipe diameter  
 Preissmann slot width  
 Domain length  
 Bump length  
 Middle of bump  
 Minimum bottom elevation  
 Maximum bottom elevation  
 Cell size  
 Initial water surface elevation  
 Prescribed total discharge at the upstream boundary   
 For subcritical flow (T3a)  
 For hydraulic jump (T3b)  
 For pressurized flow (T3c)  
 Prescribed water depth at the downstream boundary   
 For subcritical flow (T3a)  
 For hydraulic jump (T3b)  
 For pressurized flow (T3c) Free-outlet 
View Large

The simulation results of mixed pipe flow are presented in Figures 6 and 7 for pressurized flow. Numerical results show that the proposed model provides better estimations of the computed water surface elevation than the ‘Sanders & Bradford (2010)’ model, including upstream water surface elevation, shock position and transition between free-surface and pressurized flow. Regarding the water discharge, transition between free-surface and pressurized flow has caused additional oscillations for both models. For the ‘Sanders & Bradford (2010)’ model, the errors accumulated around the bump can still be observed for mixed or pressurized flow due to the adoption of HLL solver.

Figure 6
Comparison of simulated water surface elevation (left), zoomed comparison in rectangular (middle) and total discharge (right) of mixed flow over a bump.
View largeDownload slide

Comparison of simulated water surface elevation (left), zoomed comparison in rectangular (middle) and total discharge (right) of mixed flow over a bump.

Figure 6
Comparison of simulated water surface elevation (left), zoomed comparison in rectangular (middle) and total discharge (right) of mixed flow over a bump.
View largeDownload slide

Comparison of simulated water surface elevation (left), zoomed comparison in rectangular (middle) and total discharge (right) of mixed flow over a bump.

Close modal
Figure 7
Comparison of simulated water surface elevation and total discharge of pressurized flow over a bump.
View largeDownload slide

Comparison of simulated water surface elevation and total discharge of pressurized flow over a bump.

Figure 7
Comparison of simulated water surface elevation and total discharge of pressurized flow over a bump.
View largeDownload slide

Comparison of simulated water surface elevation and total discharge of pressurized flow over a bump.

Close modal

Dambreak with bottom discontinuity

Listen

This test case aims at validating the convergence ability of the model on a transient dambreak problem over a discontinuous bottom with no friction. Two cases are tested with different bottom topography: step-up case (T4a) and step-down case (T4b). The domain length is with a mesh of uniform cell ⁠. The bottom step (discontinuity) is positioned at ⁠. For the step-up case (T4a), the bottom elevations at left and right sides of the discontinuity are and ⁠, respectively. For the step-down case (T4b), the bottom elevations at left and right sides of the discontinuity are and ⁠, respectively. The dam is located at ⁠, and water velocity is zero everywhere at ⁠. For the test case T4a, the water surface levels at the left and right sides of the dam are and ⁠, respectively. For the test case T4b, the water surface levels at the left and right sides of the dam are and ⁠, respectively. The dam-break flows instantaneously start at ⁠. The computational results with different mesh sizes are compared at when the wave fronts are not reached the domain boundaries. The simulation parameters of this test are summarized in Table 3. The pipe diameter is ⁠.

Table 3

Parameters for the dambreak problem with bottom discontinuity (T4)

SymbolMeaningValue
 Pipe diameter  
 Domain length  
 Cell size  
 Bottom elevation on the left-hand side of the dam  
 Bottom elevation on the right-hand side of the dam  
 For test case T4a  
 For test case T4b  
 Initial water depth on the left-hand side of the dam  
 Initial water depth on the right-hand side of the dam  
 Computation time  
SymbolMeaningValue
 Pipe diameter  
 Domain length  
 Cell size  
 Bottom elevation on the left-hand side of the dam  
 Bottom elevation on the right-hand side of the dam  
 For test case T4a  
 For test case T4b  
 Initial water depth on the left-hand side of the dam  
 Initial water depth on the right-hand side of the dam  
 Computation time  
View Large

The simulation results at of the two test cases with bottom stepping up or down, respectively, are presented in Figure 8. It shows that, for both test cases, the solutions converge as mesh is refined. No numerical instability is produced at the position of the bed discontinuity during the dam-break process. Notice that for the step-up case T4a, a spike is observed at ⁠, this may be due to that the contact wave is ignored in the HLL Riemann solver adopted in the current study.

Figure 8
Convergence test of dambreak flow with bottom discontinuity.
View largeDownload slide

Convergence test of dambreak flow with bottom discontinuity.

Figure 8
Convergence test of dambreak flow with bottom discontinuity.
View largeDownload slide

Convergence test of dambreak flow with bottom discontinuity.

Close modal

Transient mixed pipe flow of dambreak

Listen

In this test case, the capability of the proposed method to reproduce transient flow over complex topography is tested by a mixed pipe flow experiment performed by Aureli et al. (2015), involving friction term. The experimental settings are illustrated in Figure 9. The red line represents the initial flow conditions with a pressure head of ⁠, and blue lines are the six gauges from G1 to G6 as illustrated in Table 4. The black line represents the pipe's bottom and top. The total length of the pipe is ⁠, and the first part of the pipe is with slope ⁠. The length of the second pipe is with slope ⁠. The inner diameter is ⁠, and the wall thickness is ⁠. At ⁠, the gate from the upstream inlet opens and then the water flows out. The left end of the pipe is partly closed to maintain the initial water height, and the right side is fully open to maintain a ventilated system.

Table 4

Location of the measuring gauges (⁠ is the distance along the pipe axis)

Pressure Gauges
G1 1.00 
G2 3.00 
G3 4.50 
G4 6.80 
G5 7.32 
G6 8.52 
Pressure Gauges
G1 1.00 
G2 3.00 
G3 4.50 
G4 6.80 
G5 7.32 
G6 8.52 
View Large
Figure 9
Bottom topography and initial flow conditions of Aureli et al. (2015). Please refer to the online version of this paper to see this figure in color: doi:10.2166/hydro.2021.155.
View largeDownload slide

Bottom topography and initial flow conditions of Aureli et al. (2015). Please refer to the online version of this paper to see this figure in color: doi:10.2166/hydro.2021.155.

Figure 9
Bottom topography and initial flow conditions of Aureli et al. (2015). Please refer to the online version of this paper to see this figure in color: doi:10.2166/hydro.2021.155.
View largeDownload slide

Bottom topography and initial flow conditions of Aureli et al. (2015). Please refer to the online version of this paper to see this figure in color: doi:10.2166/hydro.2021.155.

Close modal

This experiment mainly aims to study simple transient flow behavior in a single sloped Plexiglas pipe as the effect of the boundary conditions has been ignored. The numerical model proposed and the method developed by Sanders & Bradford (2010) is used to reproduce the experimental results. The cell size ⁠, and time step ⁠. Manning coefficient is set to ⁠, and wave celerity is chosen as with ⁠. Transmissive and closed boundaries are used in this work (Toro 2001). The numerical results are compared with the measurement data illustrated in Figure 10. The horizontal line in each graph is the pipe crown.

Figure 10
Comparison of mixed flow in circular pipe between experimental data and numerical results.
View largeDownload slide

Comparison of mixed flow in circular pipe between experimental data and numerical results.

Figure 10
Comparison of mixed flow in circular pipe between experimental data and numerical results.
View largeDownload slide

Comparison of mixed flow in circular pipe between experimental data and numerical results.

Close modal

As presented in Figure 10, compared with the experimental data, results show that two methods tackling slope source term can both reproduce transient mixed pipe flow of dambreak, and the curve fitting between them is satisfied. Differences can also be identified that simulated pressure is higher than the experiment in gauges G4 and G5 and lower in G6. This may be due to the minimum water depth used and adopted for free-surface flow calculation. One can observe that both the proposed scheme and Sanders & Bradford (2010) show satisfying agreement with experimental measurements.

Listen

In this paper, an augmented HLL Riemann solver is extended to circular cross-sections and mixed pipe flows over complex topography. The proposed model provides accurate predictions in reproducing steady flows with different hydraulic flow regimes. The summarized main features of this scheme are presented as follows:

  • An augmented HLL Riemann solver has been extended from rectangular shape to the currently studied pipe flow system by projecting the source term onto the Jacobian's eigenvector basis; the source term is discretized approximately as a succession of discontinuities separating regions of constant state.

  • The proposed model performs better than the reference solution obtained with the ‘Sanders & Bradford (2010)’ method when simulating steady flows over a bump, including sub-, trans-critical and hydraulic jump.

  • More accurate estimations than the reference are generated by the proposed model for mixed and pressurized pipe flow.

  • Transient dambreak pipe flow has been reproduced by the proposed model, including wetting or drying front and the calculation of the friction term.

The effectiveness and robustness of the proposed numerical scheme are verified by comparing the simulated results with numerical or measurement results. This method can also be applied to other cross-sections like trapezoid with the calculation of Riemann invariant over complex topography.

Listen

This work is funded by the National Natural Science Foundation of China (Grant No. 51922096) and the Excellent Youth Natural Science Foundation of Zhejiang Province, China (LR19E080003). The authors warmly thank Pierre-André Garambois, INRAE, for his final lecture of the manuscript.

Listen

All relevant data are included in the paper or its Supplementary Information.

Abbott
M. B.
Cunge
J. A.
1982
Engineering Applications of Computational Hydraulics: Hommage to A. Preismann. Vol. II: Numerical Models in Environmental Fluid Mechanics. Pitman Advanced Publishing Program, Boston, MA and London
.
Arrault
A.
Finaud-Guyot
P.
Archambeau
P.
Bruwier
M.
Erpicum
S.
Pirotton
M.
Dewals
B.
2016
Hydrodynamics of long-duration urban floods: experiments and numerical modelling
.
Natural Hazards and Earth System Sciences
16
(
6
),
1413
–
1429
.
Google Scholar
 
Aureli
F.
Dazzi
S.
Maranzoni
A.
Mignosa
P.
2015
Validation of single-and two-equation models for transient mixed flows: a laboratory test case
.
Journal of Hydraulic Research
53
(
4
),
440
–
451
.
Google Scholar
 
Bousso
S.
Daynou
M.
Fuamba
M.
2012
Numerical modeling of mixed flows in storm water systems: critical review of literature
.
Journal of Hydraulic Engineering
139
(
4
),
385
–
396
.
Google Scholar
 
Capart
H.
Eldho
T.
Huang
S.
Young
D.
Zech
Y.
2003
Treatment of natural geometry in finite volume river flow computations
.
Journal of Hydraulic Engineering
129
(
5
),
385
–
393
.
Google Scholar
 
Chen
S.
Garambois
P.-A.
Finaud-Guyot
P.
Dellinger
G.
Mose
R.
Terfous
A.
Ghenaim
A.
2018
Variance based sensitivity analysis of 1d and 2d hydraulic models: an experimental urban flood case
.
Environmental Modelling & Software
109
,
167
–
181
.
Google Scholar
 
Fraccarollo
L.
Toro
E. F.
1995
Experimental and numerical assessment of the shallow water model for two-dimensional dam-break type problems
.
Journal of Hydraulic Research
33
(
6
),
843
–
864
.
Google Scholar
 
Guinot
V.
2012
Wave Propagation in Fluids: Models and Numerical Techniques
.
John Wiley & Sons
,
New Jersey
.
Google Scholar
 
Harten
A.
Lax
P. D.
Leer
B. v.
1983
On upstream differencing and Godunov-type schemes for hyperbolic conservation laws
.
SIAM Review
25
(
1
),
35
–
61
.
Google Scholar
 
Hunter
N. M.
Bates
P. D.
Horritt
M. S.
Wilson
M. D.
2007
Simple spatially-distributed models for predicting flood inundation: a review
.
Geomorphology
90
(
3
),
208
–
225
.
Google Scholar
 
Kerger
F.
Archambeau
P.
Erpicum
S.
Dewals
B. J.
Pirotton
M.
2011a
An exact Riemann solver and a Godunov scheme for simulating highly transient mixed flows
.
Journal of Computational and Applied Mathematics
235
(
8
),
2030
–
2040
.
Google Scholar
 
Kerger
F.
Archambeau
P.
Erpicum
S.
Dewals
B. J.
Pirotton
M.
2011b
A fast universal solver for 1d continuous and discontinuous steady flows in rivers and pipes
.
International Journal for Numerical Methods in Fluids
66
(
1
),
38
–
48
.
Google Scholar
 
León
A. S.
Ghidaoui
M. S.
Schmidt
A. R.
García
M. H.
2006
Godunov-type solutions for transient flows in sewers
.
Journal of Hydraulic Engineering
132
(
8
),
800
–
813
.
Google Scholar
 
LeVeque
R. J.
2011
A well-balanced path-integral f-wave method for hyperbolic problems with source terms
.
Journal of Scientific Computing
48
(
1–3
),
209
–
226
.
Google Scholar
 
Li
Q.
Liang
Q.
Xia
X.
2020
A novel 1d-2d coupled model for hydrodynamic simulation of flows in drainage networks
.
Advances in Water Resources
137
,
103519
.
Google Scholar
 
Liang
Q.
Borthwick
A. G.
2009
Adaptive quadtree simulation of shallow flows with wet–dry fronts over complex topography
.
Computers & Fluids
38
(
2
),
221
–
234
.
Google Scholar
 
Liang
Q.
Marche
F.
2009
Numerical resolution of well-balanced shallow water equations with complex source terms
.
Advances in Water Resources
32
(
6
),
873
–
884
.
Google Scholar
 
Maranzoni
A.
Dazzi
S.
Aureli
F.
Mignosa
P.
2015
Extension and application of the Preissmann slot model to 2d transient mixed flows
.
Advances in Water Resources
82
,
70
–
82
.
Google Scholar
 
Martins
R.
Leandro
J.
Chen
A. S.
Djordjević
S.
2017
A comparison of three dual drainage models: shallow water vs local inertial vs diffusive wave
.
Journal of Hydroinformatics
19
(
3
),
331
–
348
.
Google Scholar
 
Mignot
E.
Paquier
A.
Haider
S.
2006
Modeling floods in a dense urban area using 2d shallow water equations
.
Journal of Hydrology
327
(
1
),
186
–
199
.
Google Scholar
 
Murillo
J.
García-Navarro
P.
2010
Weak solutions for partial differential equations with source terms: application to the shallow water equations
.
Journal of Computational Physics
229
(
11
),
4327
–
4368
.
Google Scholar
 
Roe
P. L.
1981
Approximate Riemann solvers, parameter vectors, and difference schemes
.
Journal of Computational Physics
43
(
2
),
357
–
372
.
Google Scholar
 
Sanders
B. F.
Bradford
S. F.
2010
Network implementation of the two-component pressure approach for transient flow in storm sewers
.
Journal of Hydraulic Engineering
137
(
2
),
158
–
172
.
Google Scholar
 
Sturm
T. W.
2010
Open Channel Hydraulics
.
McGraw-Hill
,
New York
.
Google Scholar
 
Teng
J.
Jakeman
A.
Vaze
J.
Croke
B.
Dutta
D.
Kim
S.
2017
Flood inundation modelling: a review of methods, recent advances and uncertainty analysis
.
Environmental Modelling & Software
90
,
201
–
216
.
Google Scholar
 
Toro
E. F.
1992
Riemann problems and the WAF method for solving the two-dimensional shallow water equations
.
Philosophical Transactions of the Royal Society of London, Series A
338
(
1649
),
43
–
68
.
Google Scholar
 
Toro
E.
2001
Shock-capturing Methods for Free-Surface Shallow Flows
.
John Wiley
,
New York
.
Google Scholar
 
Toro
E. F.
2013
Riemann Solvers and Numerical Methods for Fluid Dynamics: A Practical Introduction
.
Springer Science & Business Media
,
New York
.
Google Scholar
 
Toro
E. F.
Garcia-Navarro
P.
2007
Godunov-type methods for free-surface shallow flows: a review
.
Journal of Hydraulic Research
45
(
6
),
736
–
751
.
Google Scholar
 
Vasconcelos
J. G.
Wright
S. J.
Roe
P. L.
2006
Improved simulation of flow regime transition in sewers: two-component pressure approach
.
Journal of Hydraulic Engineering
132
(
6
),
553
–
562
.
Google Scholar
 
Vasconcelos
J. G.
Wright
S. J.
Roe
P. L.
2009
Numerical oscillations in pipe-filling bore predictions by shock-capturing models
.
Journal of Hydraulic Engineering
135
(
4
),
296
–
305
.
Google Scholar
 
Wiggert
D. C.
1972
Transient flow in free-surface, pressurized systems
.
Journal of the Hydraulics Division
98
(
1
),
11
–
27
.
Google Scholar
 
Zheng
F.
Westra
S.
Leonard
M.
2015
Opposing local precipitation extremes
.
Nature Climate Change
5
(
5
),
389
–
390
.
Google Scholar
 
© 2021 The Authors
This is an Open Access article distributed under the terms of the Creative Commons Attribution Licence (CC BY-NC-ND 4.0), which permits copying and redistribution for non-commercial purposes with no derivatives, provided the original work is properly cited (http://creativecommons.org/licenses/by-nc-nd/4.0/)