The main contribution of the paper is to incorporate pipe-wall viscoelastic and unsteady friction in the derivation of the water-hammer solutions of non-conservative hyperbolic systems with conserved quantities as variables. The system is solved using the Godunov finite volume scheme to obtain numerical solutions. This results in the appearance of a new term in the mass conservation equation of the classical governing system. This new numerical algorithm implements the Godunov approach to one-dimensional hyperbolic systems of conservation laws on a finite volume stencil. The viscoelastic pipe-wall response in the mass conservation part of the source term has been modeled using generalized Kelvin–Voigt theory. For the momentum part of the source term a fast, robust and accurate numerical scheme linked to the Lambert W-function for calculating the friction factor has been used. A case study has been used to illustrate the influence of the various formulations; a comparison between the classical solution, the numerical solution including quasi-steady friction, the numerical solution incorporating the viscoelastic effects, and measurements are presented. The inclusion of viscoelastic effects results in better agreement between the measured and solved values.
Incorporate the viscoelastic pipe and unsteady friction in water-hammer formulations.
Present the solutions of non-conservative hyperbolic systems with conserved quantities as variables.
Adapt the finite volume approach, which was developed for hyperbolic equations to tackle systems that are not in ‘proper’ conservation form.
Give a new numerical algorithm in finite volume for study several variants.
The following notation is used in this contribution.
Variables and operators
Celerity of the pressure waves
Jacobian matrix of F respecting to U
Cross-sectional pipe area
Diameter of the pipe
Pipe wall width
Young's modulus of the pipe material
Modulus of the spring of -element elasticity
Bulk modulus of elasticity
Vector flux in the x-direction
Quasi-steady part of the friction factor
Creep-compliance of the spring of Kelvin–Voigt -element
Absolute roughness of pipe
Cells number in the computational domain
Vector source term
Longitudinal fluid velocity
Unit vector in the x-direction
Calculation time step
Maximum tolerable calculation time step
Maximum tolerable calculation time step for stability of the discretized source term
Cell size in x-direction
Retarded viscoelastic strain
Constant eigenvalue of the Jacobian matrix
Mass of fluid per unit length of pipe
Retardation time of the dashpot of -element
Solution obtained at n after the solution of the conservation part in the -direction
Average value between n and
In the detailed design of pipeline systems, water-hammer analysis is very important to select pipe characteristics and to specify surge control to prevent and/or mitigate any excessive pressures. The detailed design further attempts to elaborate the complete system response using various tools such as numerical modeling. Adamkowski & Lewandowski (2006) and Shamloo et al. (2015) reviewed a quasi-steady and four unsteady friction formulations (Zielke 1968; Trikha 1975; Brunone et al. 1991; Vardy & Brown 2003) for transient pipe flow. The Trikha (1975) model gives the Zielke (1968) model with assumptions that simplify the scheme. Bousso & Fuamba (2013) have numerically represented the unsteady friction in transient two-phase flow. Other unsteady friction alternatives have been proposed by other authors (Zidouh et al. 2009; Weinerowska-Bords 2015). Their contribution bridges the gap between the detailed phase and the general concept but cannot be used as a full purpose tool for transient analysis since certain key parameters such as the effects of the viscoelastic pipe-wall are not incorporated.
Bergant et al. (2008a, 2008b) have experimentally studied key criteria that influence the pressure profile predicted by the water-hammer phenomenon. In general, numerical simulation tools assume that the friction factor is considered steady, quasi-steady, or dynamic; this has the tendency to overestimate the water-hammer peak pressures since the primary mechanism that may significantly affect pressure waveforms is the viscoelastic of the non-steel pipe-wall. In materials with viscous and elastic properties, the mechanical damping of the pressure response is larger than the damping phenomenon due to the fluid friction during the event of transient occurrences (Gally et al. 1979; Rieutord & Blanchard 1979; Rieutord 1982; Franke 1983; Güney 1983; Suo & Wylie 1990; Covas et al. 2004, 2005). Ferrante & Capponi (2018) have compared viscoelastic models with a different number of parameters for transient simulations.
The Godunov method is a conservative numerical scheme for solving PDE in CFD. This well-known conservative finite-volume approach computes a relative or exact Riemann solver at each cell. It is first-order accuracy in both time and space. This framework has been used in several papers to design first-order quasi-linear numerical methods in one-dimensional non-conservative hyperbolic systems, in the finite volume framework; see for example, Muñoz-Ruiz & Parés (2007), Guinot (2002a, 2002b, 2004), Bousso & Fuamba (2013), and Seck et al. (2017). Guinot (2003) focused on the scheme of the Godunov method to solve the classical water-hammer equations using finite volumes with steady friction. Pal et al. (2020) compared the first- and second-order FVMs based on the Godunov scheme with the method of characteristics (MOCs) and experimental water-hammer measurements available in the literature. Seck et al. (2017) developed the derivation of water-hammer formulations’ hyperbolic systems of conservation laws with unsteady friction, which are not applicable for viscoelastic pipes because their non-negligible effects have not been taken into account.
Derived from the principles of conservation of momentum and mass, the water-hammer equations integrating unsteady friction and viscoelastic pipe-wall are expressed in hyperbolic systems of conservation laws with conserved quantities as variables. The presence of the viscoelastic pipe-wall terms preclude presentation in the integral (conservation) form and cause a source term to appear. The hyperbolic systems admit weak solutions in the form of discontinuities or ‘shocks’. The treatment of these discontinuities requires special treatment of the flux function to avoid spurious oscillations.
The innovative components of this contribution are to present the water-hammer derivation equations with the effects of viscoelastic pipe-wall and unsteady friction, and adapt the finite volume approach with conserved variables, for hyperbolic equations of conservation laws. The aim of this contribution is to provide a new user-friendly complete numerical algorithm in finite volume for studying several variants or design water-hammer, in a short period, of the series of different characteristics of viscoelastic pipes in the distribution network that can be of benefit to engineers.
This paper focuses on implementation of the Godunov approach to the hyperbolic systems in one-dimensional conservation laws that describe the phenomenon of water-hammer integrating the viscoelastic pipe-wall and unsteady friction. The solution of a hyperbolic system of conservation laws by Godunov-type algorithms comprises six steps:
discretization of space in finite volumes,
construction of generalised Riemann problems (GRP) at the cell,
conversion of the GRP into equivalent Riemann problems (ERP),
solution of the ERP,
balance over the cells,
incorporation of the effect of the source terms.
The computation of the fluxes is detailed followed by a description of the complete numerical algorithm. A case study is examined with comparisons being made between the numerical solution with quasi-steady friction, the numerical solution with viscoelastic pipe-wall, a classical solution, and the experimental results.
Derivation from the conservation of mass and momentum
, is the flow conserved variable vector, the flux vector, the source term vector, D the inner pipe diameter (m), the Vardy's shear decay constant related to the flow regime, x the unit vector in x-direction, t the time (s), the mass of liquid per unit length of pipe (kg/m), the mass discharge (kg/s), the fluid density (kg/m3), the cross-sectional pipe area (m2), V the longitudinal fluid velocity (m/s), wave speed (m/s), the steady state friction coefficient, and p the pressure (Pa) where the formula will be given later (Equation (21)).
In the definitions and basic notion of conservation laws, the source term can be a function of time t, space coordinate x and flow variable vector U, combined or separately. In this form (Equation (1)), where only the conservation part remains on the left, the term on the right is considered as being part of the source term. This source term thus contains an additional term . is referred to as a hyperbolic source term, if it incorporates the derivative . Systems that can be written under the form of Equation (1) are also called systems of conservation laws with source term, or hyperbolic systems of conservation laws, or balance laws (Parés 2006). The incorporation of the unsteady friction results in the appearance of the quasilinear term in the hyperbolic system (Equation (1)) thus destroying the classical integral form of the equations (Seck et al. 2018).
The viscoelastic pipe-wall inclusion
From a physical point of view, the third term in Equation (16) is the flux in the radial direction. That component is a function of conserved vector variable , and time derivative of retarded strain , combined, but does not contain functions of any of derivatives of .
The next section describes the solution process for the system.
The first-order time splitting is used. The viscoelastic model (based on the generalized Kelvin–Voigt application) in the mass conservation part of the source term is solved using the modified Covas et al. (2005) numerical scheme. For the momentum part of the source term, the Colebrook–White equation is solved using a robust, fast, and accurate numerical scheme due to Clamond (2009).
Solutions at the internal cells and at the boundaries
Discretization of the source term
Computational time step
The user-friendly proposed algorithm is presented in Figure 2.
RESULTS AND DISCUSSION
A case study has been simulated. Data from the experimental test conducted by Covas et al. (2004) have been compared with the numerical results.
The results of experimental studies conducted by Covas et al. (2004) obtained from a PE pipe-rig at Imperial College and considered to be of good quality are compared with the mathematical model. The experimental setup is composed of a 277 m long PE pipe with a thickness of 6.25 mm and an inner diameter of 50.6 mm that connects to an upstream air vessel. Initial steady flow (Q0 = 1.008 l/s) are followed by a transient occurrence initiated by the sudden closure of a globe valve at the downstream end. = 1,000 kg/m3 and g = 9.81 m/s2. The initial wave celerity a0 was evaluated at 395 m/s. The pressure waves were computed and compared with the corresponding experimental values at location #1 (distance x= 273 m from upstream), location #5 (distance x= 117.5 m from upstream) and location #8 (distance x= 199 m from upstream). The values of the creep coefficients for the retardation times are specified as follows. A full description of the collected data and the experimental system may be found in Covas et al. (2004).
J0 = 0.7.10−9 Pa−1
J1 = 0.1394.10−9 Pa−1; = 0.05 s;
J2 = 0.0062.10−9 Pa−1; = 0.5 s;
J3 = 0.1148.10−9Pa−1; = 1.5 s;
J4 = 0.3425.10−9Pa−1; = 5 s;
Results and discussion
The water-hammer equations incorporating viscoelastic pipe-wall were solved using the Godunov scheme detailed above. In keeping with the objectives of this paper, results were obtained for the analytical solution, the numerical solution with quasi-steady friction, and the numerical solution including the tuned viscoelastic model using a time step of . The experimental measurements are compared with these results. Good accord between measured and computed viscoelastic profiles has been obtained (Figures 3–5).
Comparison between the results indicates that the inclusion of viscoelasticity has a serious impact on the pressure profiles. The quasi-steady friction model does not reproduce the progression of the aspect of the pressure waves very well. They also display phase shift errors. The quasi-steady friction model undervalues the propagation and damping predicted by the (physically) more accurate viscoelastic model which tends to attenuate the magnitude of the overpressures more rapidly. Due to damping, the waves will be of decreasing amplitude until the final equilibrium pressure is reached. In the present context given the simple geometry and non-varying fluid properties (the fluid is slightly compressible), the basic Godunov scheme is appropriate for the analysis of water-hammer phenomenon. The wave profile of the viscoelastic line is smooth and the wave fronts are smeared which justify the presence of a rarefaction wave. That front smearing is caused by the exponential term contained in the mass conservation equation (Equation (40)). The Godunov scheme introduces a non-negligible diffusion. This explains the use of a large number of calculation cells to obtain a quality solution similar to that of a higher order scheme that requires few computing cells. More cells and therefore more computational effort is required for one-order scheme. But, for a robust and fast one-order scheme that runs in a short time with many cells, the use of a higher order scheme is not justified, because it will only improve the solution just weakly.
The water-hammer equations, in conservation form, combining viscoelastic pipe-wall and unsteady friction, are developed and numerically solved using a finite volume formulation for engineering applications. A robust, simple, accurate, and fast computational algorithm based on the Godunov scheme for one-dimensional hyperbolic systems is presented in some detail. Incorporation of the viscoelastic pipe-wall results in the appearance of a new term in the mass conservation part of the source term in the hyperbolic system of governing partial differential equations. From a physical point of view, this source term in mass conservation part is the flux in the radial direction. A case study has been presented to compare the separate impacts on wave attenuation with time. The findings indicate, as expected, inclusion of viscoelasticity reduces the peak water-hammer pressures. The model developed appears to be able to satisfactorily predict transient pressures in viscoelastic pipes and can be a useful tool for qualitative analysis.
These equations, as expected, admit ‘weak’ solutions with hydraulic jumps being the discontinuities. Using local first-order approximations to handle the spurious behavior near shocks does work without the expense of additional computational time. Computational time for the case studied took only 1 minute and 28 secs on a PC. Thus, for engineering applications, several variants for the design of series of different characteristic viscoelastic pipes in distribution networks may be studied in a short period of time with a spatial discretization not necessarily uniform.
The algorithm introduced in this paper could be used however, of course with some modifications, in the high-order schemes. This first-order proposed algorithm, in the context of hydraulic transients in pressurised pipelines incorporating viscoelastic pipe-wall and unsteady friction, is an interesting start to achieve this aim.
DATA AVAILABILITY STATEMENT
All relevant data are included in the paper or its Supplementary Information.