Abstract

The numerical and analytical models used for transient simulations, and hence for the pressurized pipe system diagnosis, require the definition of a rheological component related to the pipe material. The introduction and the following widespread use of polymeric material pipes, characterized by a viscoelastic behavior, increased the complexity and the number of parameters involved in this component with respect to metallic materials. Furthermore, since tests on specimens are not reliable, a calibration procedure based on transient test is required to estimate the viscoelastic parameters. In this paper, the trade-off between viscoelastic component accuracy and simplicity is explored, based on the Akaike criterion. Several aspects of the calibration procedure are also examined, such as the use of a frequency domain numerical model and of different standard optimization algorithms. The procedure is tested on synthetic data and then it is applied to experimental data, acquired during transients on a high density polyethylene pipe. The results show that the best model among those used for the considered system implements the series of a spring with three Kelvin–Voigt elements.

INTRODUCTION

The numerical models for the simulation of transients are an important tool for the design and the diagnosis of pressurized pipe systems (Covas et al. 2004; Duan et al. 2010b, 2012; Meniconi et al. 2011; Gong et al. 2013; Massari et al. 2014; Duan 2017; Laucelli et al. 2016; Sun et al. 2016). Such models require the definition of a functional relationship between stresses and strains in the pipe material, i.e., its rheological model. The spreading of pipes of polymeric materials, such as high density polyethylene (HDPE) and polyvinyl chloride (PVC), and their use in pressurized pipe systems instead of metallic pipes, increased the complexity of the rheological models, from elastic to viscoelastic. The increase in the complexity of the models corresponded to an increase in the number of the parameters they depend on and to the need of reliable and efficient algorithms to calibrate them.

For the specific purpose of transient modeling in polymeric pipes, the evaluation on specimens of the rheological parameters, i.e., by creep and relaxation tests in laboratory, is not successful (Ghilardi & Paoletti 1986; Covas et al. 2005). For this reason, the most common procedure used so far is the calibration of the parameters of the implemented rheological models using transient tests.

In the literature, several case studies and calibration procedures have been presented. As an example, Pezzinga & Scandura (1995) used a trial and error procedure to calibrate a rheological model consisting of one Kelvin–Voigt (KV) element with two parameters, by means of transient tests on a HDPE pipe. In the following papers, Pezzinga (2000, 2002) and Pezzinga et al. (2014) calibrated the same rheological model by means of a micro-genetic algorithm. The calibration procedure is based on the maximization of a function evaluated as the inverse of the sum of the squared residuals (RSS) between simulated and experimental data. The numerical model is obtained by the integration of the governing equations with the method of the characteristics (MOC).

Covas et al. (2005) used the experimental data collected from a PE pipe-rig to calibrate a rheological model composed by a different number of KV elements, fixing some of the parameters and varying the others. The calibration is based on the minimization of the RSS by means of a Levenberg–Marquardt algorithm. A sensitivity analysis on the RSS suggests that the optimal results are obtained with five KV elements and five calibration parameters. The need of a proper choice of the time step in the MOC model is discussed.

Soares et al. (2011) used a genetic algorithm to calibrate different rheological models on data coming from tests carried out in an experimental hydraulic circuit composed of PVC pipes. Numerical MOC models implementing different rheological models, with up to three KV elements, are calibrated minimizing the RSS by means of both a genetic and a Levenberg–Marquardt algorithm. The algorithms are used in series since the genetic algorithm is considered as a good tool for the detection of local minima but not accurate enough for the global minimum determination. The considered optimal value is obtained by one KV element and fixing some parameters.

Although in the available literature different choices are made with regard to the complexity of the model, the number of the involved parameter is not considered as a parameter itself. In other words, the complexity level of the used viscoelastic models is not investigated as a variable of the problem and determined by means of an objective and rational criterion.

In this paper, several aspects of the calibration of the viscoelastic parameters are investigated. In the first part, the implementation of the viscoelastic elements in a frequency domain model is described. The choice of the optimization function and algorithm is discussed, giving some general remarks regarding the calibration procedure. A method for the estimation of the appropriate number of the viscoelastic model parameters is then presented. The proposed approach is applied to synthetic data, obtained corrupting with a white noise the numerical results of an implemented viscoelastic model, and then to experimental data, obtained by tests carried out on a HDPE pipe. The obtained results and the original contributions of the paper are summarized in the conclusions.

VISCOELASTIC RHEOLOGICAL MODEL

The rheology of a pipe material can be described introducing simple components, corresponding to elementary behaviors, and combining them to describe more complex behaviors. A linear viscoelastic material, being in between a linear elastic and a linear viscous material, can be modeled as a combination of elementary linear elastic (springs) and linear viscous (dashpots) elements. The combination of a spring and a dashpot in parallel is the so-called KV element. It can be analytically described as:  
formula
(1)
where the first and second terms on the right represent the spring and the dashpot, respectively. In Equation (1), are the stresses and are the strains, E denotes the Young's module of the spring and the viscosity of the dashpot, t is the time while the subscript denotes the model.
A KV element coupled in series with a spring is usually referred to as standard linear solid (S1). More complex systems can be obtained putting in series with the spring not just a single KV element but a series of m KV elements (Sm). In this case, it is:  
formula
(2)
where and and the subscript 0 denotes the spring in series with the KV elements.

Many authors (e.g., Ghilardi & Paoletti 1986, 1987; Pezzinga & Scandura 1995; Covas et al. 2005; Lee et al. 2009; Duan et al. 2012; Keramat et al. 2012; Pezzinga 2014; Evangelista et al. 2015) have presented the implementation of S rheological components in MOC numerical models, to relate hoop stresses and strains in the pipe. Other authors have presented the implementation of the same rheological components in frequency domain integration models (e.g., Duan et al. 2010a, 2012; Vítkovský et al. 2011; Lee et al. 2013; Capponi et al. 2017a; Ferrante & Capponi 2017b). With respect to the time domain, the frequency domain integration is faster and easy when the variation in time at one section of the system is sought. The frequency domain integration requires a linearization of the governing equations but limited to the second order of the steady-friction resistance term, when the boundary conditions are linear (Capponi et al. 2017c). Since the viscoelastic term in the governing equations can be represented as a convolution in time (Weinerowska-Bords 2015) it does not require a linearization and, due to the well-known properties of Fourier and Laplace transforms of such integrals, it is easier to be handled in the frequency domain than in the time domain. Similarly, the unsteady-friction term is a convolution integral in the time domain as well (Weinerowska-Bords 2015) and can be considered and easily implemented in the frequency domain models. With reference to the considered conditions and to the results presented in the following, the effects of viscoelasticity prevail on those introduced by the friction term (Duan et al. 2010a, 2012). Hence, for the sake of simplicity and considering that the main purpose of this paper is the comparison of the viscoelastic models, a frequency domain numerical model without an unsteady-friction term is used.

THE UNSTEADY-FLOW MODEL IN THE FREQUENCY DOMAIN

For the pipe material rheological model S, defined by Equation (2), the total hoop strain is the sum of the spring strain, and of the strain associated with the remaining part of the model, . Due to the series arrangement, the stress applied to the spring is the same applied to the remaining part of the model. Under this and other common assumptions (Wylie & Streeter 1993; Chaudhry 2014), the Fourier transform of the linearized equations governing the one-dimensional transient flow can be obtained:  
formula
 
formula
(3)
where and is the angular frequency. The dependent variables, h, q, and , are the Fourier transform of the pressure head, H, flow, Q, and pipe hoop strain, , respectively. In Equation (3), the capacitance, , the inertance, , and the resistance, , are introduced, where is the friction factor, g is the gravitational acceleration, A and D are the pipe cross-sectional area and diameter. Only the strain component of the spring, , contributes to the evaluation of the wave speed, a, while the term in Equation (3) takes into account the characteristics of the remaining part of the rheological model, i.e., the m introduced KV elements. Combining Equation (3) with the equations of the rheological model, two symmetric equations with the dependent variables h and q can be derived by simple algebraic manipulations:  
formula
 
formula
(4)
where  
formula
 
formula
(5)

The function depends on the chosen rheological model of the pipe material (Ferrante & Capponi 2017b). In the fifth column of Table 1 the formulae are given for the S models, for different values of m. In these formulae, assuming that the thickness of the pipe wall , and , it is , where is the pipe constraint coefficient and is the water density. In the same table the model parameters are also specified. Since all the used models require a spring in series with the rest of the element arrangement, the Young modulus of this spring, , should be added directly to the list or indirectly, i.e., considering the wave speed a as a further parameter.

Table 1

Used viscoelastic models

Model Sketch  Parameters f(ω
EL  E0 
S1  E0, ES1, ηS1  
S2  E0, ES1, ηS1,ES2, ηS2  
… … … … … 
Sm  2m + 1 E0, ES1, ηS1,…, … ESm, ηSm  
Model Sketch  Parameters f(ω
EL  E0 
S1  E0, ES1, ηS1  
S2  E0, ES1, ηS1,ES2, ηS2  
… … … … … 
Sm  2m + 1 E0, ES1, ηS1,…, … ESm, ηSm  
The integration of Equation (4) for the case of a reservoir-pipe-valve (RPV) system yields the expression (Wylie & Streeter 1993; Ferrante & Brunone 2003):  
formula
(6)
where is the impedance function at the pipe downstream end, that is just upstream of the valve, , and is the length of the pipe.

In Equation (6), is introduced as a parameter, without any need of further space discretization and compatibility conditions between space and time grids. If a reduced number of measurement sections is used in the simulation, the large amount of calculations needed by the time domain integration, at each time step and for each node of the grid, can be strongly reduced by the frequency domain integration. Since for the calibration procedures the dependent variables are evaluated at a few sections, frequency domain models are much more efficient than their time domain correspondent. The reliability and the accuracy of frequency domain solutions rely on other parameters (Lee & Vítkovský 2010; Lee 2013).

With respect to the approximations introduced by the steady-friction term linearization, Capponi et al. (2017c) have shown that, with reference to the test conditions of this investigation, their effects are negligible when compared to those due to viscoelasticity. Furthermore, the computational burden is reduced of about two orders of magnitude, with respect to the time domain (e.g., Capponi et al. 2017a, 2017b; Ferrante & Capponi 2017b). Hence, considering that the increase in the simulation efficiency justifies the introduced approximations, the frequency domain model is used in the following.

MODEL CALIBRATION

The values of the parameters of the numerical model can be evaluated by means of a calibration procedure, i.e., by the comparison of the results of the numerical model with the experimental data. For a given model, the values of the parameters are varied to minimize the distance between the model and the actual phenomenon that produced the experimental data. Three different aspects can be considered: (i) the definition of a measure for the distance between the chosen model and the actual phenomenon, (ii) the algorithm to minimize such a distance, and (iii) the definition of a measure to compare the distance from the actual phenomenon of models with a different number of parameters.

The measure of the distance between a model and the actual phenomenon

Several functions can be used to measure the distance between the actual phenomenon and the model, or at least the distance between the acquired pressure head values during a transient, and the corresponding simulated values, H, produced by a model (Savic et al. 2009). The calibration procedures aim to estimate the parameter values that minimize or maximize the value of this function, defined as the optimization function.

For transients in pressurized pipe systems, the sum of squares of the residuals, RSS, can be used as the optimization function (e.g., Covas et al. 2005; Jung & Karney 2008; Savic et al. 2009):  
formula
(7)
where denotes the sample order and n is the total number of samples. Other functions based on RSS can also be used, such as the inverse (Pezzinga 2014; Pezzinga et al. 2014):  
formula
(8)
or a weighted RSS (Kapelan et al. 2003):  
formula
(9)
where are the weights.
The Nash & Sutcliffe (1970) coefficient, NS:  
formula
(10)
where is the mean value of , is the most used efficiency index in hydrology, when measured and forecasted data are compared and it has been also used for transient tests in pressurized pipes (Ferrante et al. 2015; Ferrante & Capponi 2017a).
In this study, an optimization function based on RSS is used:  
formula
(11)
because it has some interesting properties for the comparison of models with a different number of parameters.

The optimization algorithm and the comparison of calibrated models with a different number of parameters

Various algorithms can be used to minimize the optimization function in the parameter space. The optimization by trial and error is not a reliable and reproducible process while the optimization by brute force, i.e., by the direct scrutiny of the function on a grid of the parameter space, is numerically cumbersome and can be used only for efficient models with a reduced number of parameters. A genetic algorithm and two other nonlinear minimization techniques, i.e., the Nelder–Mead or the Levenberg–Marquardt algorithms were used for the presented calibration procedure. All the mentioned algorithms are available as libraries of interpreters and compilers. A freeware collection of minimization functions is provided in Octave (Eaton et al. 2015).

The minimization of the optimization function leads the chosen model as close as possible to the actual phenomenon but it does not provide enough information about the relative distance of different models from the actual phenomenon, at least when the comparison is among calibrated models with a different number of parameters. In this case, the measure of the relative distance cannot be solely based on the residual analysis and on RSS but needs also to consider the number of the parameters involved.

The Akaike criterion is a rigorous model selection criterion based on the Kullback–Leibler information (Burnham & Anderson 2003). It provides a measure of the relative distance between models and the unknown actual phenomenon that takes into account the number of the parameters. In the case of least square estimation with normally distributed residuals it can be written as:  
formula
(12)
where n is the total number of samples. The trade-off between the minimization of and the minimization of the number of the parameters, , that is between underfitting and overfitting, is expressed by the balance of the two terms of Equation (12). The measure provided by the is relative, in the sense that it can only be used to compare different models in a chosen set of models, fitting the same data. It does not provide an absolute value for just one model and cannot be used to compare models fitted to different experimental data.

The model assessed to be the best among the others in a set is the one with the minimum value of , i.e., . For all the other fitted models, the larger the value of , the less plausible it is that they can be considered as the best model. As a rule of thumb (Burnham & Anderson 2003), all the models with have still substantial support to be estimated as the best, models with have considerably less support while models with have no essential support.

It is worth noting that Equation (12) depends on the value of . If the parameter estimation procedure is carried out perfectly, the obtained value of for a chosen model is the minimum. Hence, the efficiency and reliability of the calibration algorithm in reaching the minimum of the objective function plays a crucial role in the comparison of the models. This aspect is of great importance for the fitting of the viscoelastic models in the case of diagnosis by means of transients, because of the low gradients of the optimization functions in the parameter space region close to the minimum. This flat region implies that the sensitivities of the optimization function with respect to the parameters are low, and hence, that the minimization algorithm has to be efficient since large variations in the parameters correspond to small variations in the optimization function values.

SYNTHETIC DATA

To verify the reliability of the calibration procedure and of the proposed criterion of model comparison, a synthetic data set was generated by means of a numerical model based on Equation (6). Two RPV systems (Figure 1), with pipes of a total length, , of 102.58 m, and diameter, D of 96.8 mm, were considered, made of two different pipe materials, similar to HDPE. For the first material, A, the rheology was considered completely defined by an S3 model with = 353.73 m/s, = 4.500 109 Pa, = 4.500 1010 Pa, = 4.500 1011 Pa, and = 9.000 109 Pa s. As a result, the KV elements were characterized by = 0.020 s, = 0.200 s, and = 2.000 s. The second material, B, was still defined by an S3 model and had the same value of a as the first one, but = 4.500 1010 and = 9.00 108 Pa s, = 9.00 109 Pa s, and = 9.000 1010 Pa s. As a result, the KV elements were characterized by the same values of of the previous model, i.e., = 0.020 s, = 0.200 s, and = 2.000 s.

The transients in the two systems were produced by the closure maneuver of the downstream valve DV starting from the initial steady-state conditions with the discharge, = 3.64 l/s. The flow variation in time at the valve, , was assumed to follow the sigmoid function:  
formula
(13)
with = 32 s and = −5.5.
Figure 1

The considered set-up: R = air vessel, PT = pressure transducer, SG = strain gage, FM = flow meter, MV = remotely controlled butterfly valve, DV = hand operated ball valve.

Figure 1

The considered set-up: R = air vessel, PT = pressure transducer, SG = strain gage, FM = flow meter, MV = remotely controlled butterfly valve, DV = hand operated ball valve.

The steady-friction term was evaluated during the simulation with the Blasius formula for turbulent flow in smooth pipes. The pressure head variations in time, or pressure signals, generated in the two systems with the two different pipe materials, were then corrupted with a white noise. The same generated noise was used for both signals, with a zero mean and a standard deviation of 10−1 m, to simulate the acquisition by an actual pressure transducer. The sample variance of the generated random noise was 1.005 10−2 m2.

Assuming that the generated synthetic signal, affected by the superimposed random noise, represents the truth and that the model based on the implemented S3 contains all and only the actual governing equations of the phenomenon, a blind calibration procedure was applied as in cases where the actual model is unknown. Since, in reality, also the number of the parameters is considered as unknown, models with a different number of KV elements and parameters were considered.

The calibration procedure for the five models (S1, S2, S3, S4, and S5) and the two pipe materials (A and B), required the use of a nonlinear algorithm to minimize . The first 20 seconds of the pressure signal were considered for the evaluation of . The upper and lower bounds of 1012 and 106 were used for the viscoelastic parameters. A genetic algorithm, with populations of 200 individuals, a number of populations equal to the number of the parameters and up to 100 generations, produced values of always greater than 3.0, far from the known solution. As a consequence, the heuristic Nelder–Mead algorithm was used instead, with a minimum tolerance of 10−16 on the function and on the parameter variations and a maximum number of 1,000 function evaluations. Since references to the use of the Levenberg–Marquardt algorithm can be found in the literature (Covas et al. 2005; Soares et al. 2011), a check was made to verify that the Nelder–Mead solution was the same as the Levenberg–Marquardt algorithm, with the same tolerances. At the end of the calibration procedure, the optimal set of parameters and the minimized value of were available for each combination of model and material.

In Figure 2 the synthetic signals are compared to those simulated implementing the viscoelastic models S1, S2, S3, S4, and S5 and considering the system made of the two materials, A and B (Figure 2(a) and 2(b), respectively). Only small differences can be appreciated between the synthetic signals and those simulated by S1 and S2, while other simulated signals overlap for both materials. Differences among results of the models S1, S2, and S3 are much more evident in Figure 3, where the absolute values of the functions of Equation (5) are plotted for the pipe materials A and B (Figure 3(a) and 3(b), respectively). The results of the models S4 and S5 were very close to those of S3 and are not plotted.

Figure 2

The synthetic pressure signal at the valve (Synth.) for the pipe materials A (a) and B (b) are compared with the results of the fitted models S1, S2, and S3. Results of S4 and S5 cannot be distinguished from the synthetic signal and are not plotted.

Figure 2

The synthetic pressure signal at the valve (Synth.) for the pipe materials A (a) and B (b) are compared with the results of the fitted models S1, S2, and S3. Results of S4 and S5 cannot be distinguished from the synthetic signal and are not plotted.

Figure 3

The function for the synthetic data and estimated by the fitting of the numerical model for the pipe materials A (a) and B (b).

Figure 3

The function for the synthetic data and estimated by the fitting of the numerical model for the pipe materials A (a) and B (b).

The wave speed was included as a parameter in the calibration procedure. The calibrated values of a were between 353.2 and 354.1 m/s, and hence with differences of less than ±0.5 m/s with respect to the value, used to generate the synthetic data.

The minima reached for decrease with the number of KV elements, m, and tend towards the sample variance of the added random noise.

Figure 4 presents the calibration results in the parameter space. Due to their significance, J and are chosen instead of E and . Dashed lines and solid circles denote the true solutions.

Figure 4

The values of and J estimated by the calibration of the numerical models for the pipe materials A (a) and B (b).

Figure 4

The values of and J estimated by the calibration of the numerical models for the pipe materials A (a) and B (b).

With reference to the values, we can conclude that, for both materials, values of S1 and S2 fall in the range of the variation of the true values. For the model S3, the values corresponding to the obtained minimum are very close to the true ones in the case of material A (Figure 4(a)) while for material B both and are clearly far from the true ones (Figure 4(b)). For S4 and S5, the missing solutions overlap with others or are outside the limits of logarithmic axis.

No correlation can be clearly indicated between the calibrated values of and the maneuver duration (about 0.15 s), the pipe characteristic time (0.58 s), and the duration of the synthetic signals used for the calibration (20 s).

Figure 4(a) shows an interesting correlation between estimated and J values, which leads the calibration to end along the line defined by the true solutions in the parameter space.

The variation of with the number m of the KV elements of the used models is shown in Figure 5. While is monotonically decreasing with m, the minima of are obtained for , corresponding to the value of the true model for both materials, A and B. The models S1 and S2 are characterized by a and hence can reasonably be excluded if compared to S3. The obtained values of , close to 4 and 7 for models S4 and S5, respectively, suggest that these models are significantly more distant from the true model than S3.

Figure 5

The variation of with the number m of the KV elements of the model used for the calibration. The crosses denote the theoretical values obtained by = 1.005 10−2 m2. For it is .

Figure 5

The variation of with the number m of the KV elements of the model used for the calibration. The crosses denote the theoretical values obtained by = 1.005 10−2 m2. For it is .

The crosses in Figure 5 denote the theoretical values obtained by assuming = 1.005 10−2 m in Equation (12), i.e., considering that the differences between calibrated and actual data are completely explained by the random noise. A calibration that does not reach the theoretical minimum produces a higher value of and affects the measure of the distance from the true model. In this sense, a reliable calibration procedure is crucial. Since hollow circles and crosses coincide in Figure 5, the used algorithms, with the given tolerances and settings, can be considered sufficiently accurate for the aims of the comparison.

To reduce the number of the parameters, some authors (Covas et al. 2005; Soares et al. 2011) prefer to fix the values of and calibrate the others. To assess the use of this strategy, five reduced S models, with fixed values of , were also considered in the set of model candidates, with m ranging from 1 to 5. The fixed values of were the first m values of 0.01, 0.10, 1.00, 10.00, and 100.0 s. In this case, the obtained value of for the reduced S2 model, with fixed values of = 0.01 s and = 0.1 s, was less than that obtained by the best calibration of the S3 model, estimating all the parameters. This result suggests that, following this strategy, the number of parameters could be further decreased.

EXPERIMENTAL DATA

The calibration procedure based on the Nelder–Mead algorithm and the model comparison methodology verified on the synthetic data were applied to the experimental data to determine the optimal number of the KV elements needed to model transients in an actual HDPE pipe.

The experimental apparatus

The tests were carried out at the Water Engineering Laboratory of the University of Perugia, Italy, on the same RPV system used to generate the synthetic data, with the same and D. The pipe was a HDPE DN110 PN10, according to UNI EN 12201 and UNI EN ISO 15494. An electromagnetic flowmeter was used to measure the discharge in the steady-state initial conditions, with an accuracy of 0.2% of the measured value. A piezo resistive pressure transducer, PT, with a full scale of 6 bar and an accuracy of 0.15% of the full scale was used to measure the pressure upstream of the maneuver valve.

At the downstream end of the pipe, a hand-operated ball valve, DV, discharging into the air and a remotely controlled butterfly valve, MV, immediately upstream of DV (Figure 1) were used to have an initial condition of = 3.64 l/s and to produce the maneuver of Equation (13), with the same parameters used for the generation of the synthetic data. The obtained variation in time of the pressure head at PT, , is shown in Figure 6.

Figure 6

The experimental pressure signal at the valve (Exp.) is compared with the numerical signals obtained by the calibrated models S1, S2, and S3. The first five seconds of (b) are shown in (a).

Figure 6

The experimental pressure signal at the valve (Exp.) is compared with the numerical signals obtained by the calibrated models S1, S2, and S3. The first five seconds of (b) are shown in (a).

The calibration of the models by the experimental data

Six models differing for the number m of the KV elements were implemented and calibrated to the experimental data. The Nelder–Mead algorithm was used to evaluate the parameters of the models minimizing the value of , with the same tolerances and limits used for the synthetic data. The comparison of the simulated and experimental data in Figure 6 and in Figure 7 confirms the results obtained for the calibrations by synthetic data. While the differences in the pressure signal are small, differences in the function are evident, at least for models S1, S2, and S3. Results of the S4, S5, and S6 models, not shown in the figures, are very close to the results of the S3 model.

Figure 7

The function obtained by the calibration of the S1, S2, and S3 models on the experimental data of Figure 6.

Figure 7

The function obtained by the calibration of the S1, S2, and S3 models on the experimental data of Figure 6.

The values of the parameters obtained for the six models by the calibration procedure are shown in Figure 8.

Figure 8

The estimated value of and obtained by the calibration of the S models on the experimental data of Figure 6, with the number of KV elements, m, ranging from 1 to 6.

Figure 8

The estimated value of and obtained by the calibration of the S models on the experimental data of Figure 6, with the number of KV elements, m, ranging from 1 to 6.

As expected, the values of obtained at the end of the calibration for each model decrease monotonically with the number of used KV elements.

The values of shown in Figure 9 allow assessment of whether the decrease is justified by the increase in the model complexity. Based on , S3 is the model closest to the actual phenomenon that generated the experimental data, among the set of models containing S1, S2, S3, S4, S5, and S6. No other model in the set shows a and hence there is no substantial support to estimate another model in the set as the best one.

Figure 9

The variation of with the number m of the KV elements of the calibrated models. For it is .

Figure 9

The variation of with the number m of the KV elements of the calibrated models. For it is .

Assuming that the value of of the model S3 corresponds to that of the noise of the experimental data, the values of can be estimated for S4, S5, and S6. The comparison of these values (crosses) with the data obtained by the calibration confirms that the increase of the number of the parameters, for , cannot be justified by the decrease in . On the other hand, the decrease in from S1 and S2 to S3 justifies the increase in the number of the parameters.

DISCUSSION OF THE RESULTS

The results shown allow some remarks on the calibration of viscoelastic models for transient test in polymeric pipes. Several issues need to be considered to achieve reliable results.

The first issue is related to the choice of the viscoelastic model. Although other models could be considered (Ferrante & Capponi 2017b), the most common choice is a series of one or more KV elements with a spring. A set of viscoelastic models with up to six KV elements is considered in this paper.

The parameters of the models can be estimated by means of transient tests, since creep and relaxation tests on specimens do not provide reliable results. This leads to the second issue, that is to define a measure of the estimation reliability of a set of model parameters. This measure is usually expressed as a function of the differences between simulated and experimental data. The estimate of the parameters requires the minimization of this function. The sum of the squares of the residuals is one of the most common choices and is used in this paper. Unfortunately, this function presents wide flat regions in the parameter space close to the optimal solutions and this leads to the third issue, that is the choice of the algorithm to be used to estimate the parameters minimizing the distance between simulated and measured data. The most used algorithms in the literature are also applied in this paper to synthetic and experimental data. The genetic algorithm with the shown settings is confirmed to be not reliable and efficient in reaching the global optimal solution. In contrast, the Nelder–Mead algorithm is able to attain the optimal solution in a reasonable computational time, even for low tolerances. This algorithm is simpler and more robust than the Levenberg–Marquardt algorithm, at least for the considered cases.

All the models in a set, with a different number of viscoelastic elements, provide optimal solutions and this leads to the fourth issue, that is, the comparison of models with a different number of parameters. The need of the trade-off between the reduction of the squared sum of the residuals and the increase in the number of the parameters can be handled using different approaches. The Akaike information criterion is used in this paper to introduce into the calibration procedure the number of the parameters as a further parameter, although other techniques can also be used for an efficient investigation of the behavioral model space (Beaumont et al. 2002; Sadegh & Vrugt 2014). The use of this criterion is tested on synthetic data and is applied to estimate the best model for experimental data. The obtained results show that the criterion can be successful and that a viscoelastic model with three KV elements, S3, is the best estimated model in the set of models with up to six KV elements.

The calibrations presented in this paper for the experimental data did not show any correlation of with the duration of the maneuver, the duration of the acquisition, and the characteristic time of the system. For the synthetic data, these durations are taken into account in other model components and are not correlated to the parameters in the used viscoelastic model equations. The difficulties in the calibration procedure suggest that the sensitivity of the sum of the squared residuals to is low. Variations of one order of magnitude of can cause small variations in the sum of the squared residuals. Probably due to this low sensitivity, the strategy of a reduction of the number of parameters by fixing values of can be successful, as for the shown application of a reduced S2 model. The reduction of the number of the parameters by fixing the values with a different order of magnitude (e.g., 0.01, 0.10, 1.00, 10.00 s) is based on a subjective choice of the fixed values and unfortunately it does not allow a general conclusion to be drawn here.

With reference to the applications, different model uses, such as for the preliminary design of a simple pipe or for the diagnosis of a complex system, require different levels of accuracy. In this paper, the obtained minima are comparable to the accuracy of a pressure transducer, as required by the diagnosis of pressurized pipe systems by transient tests. Furthermore, the initial 20 s of the experimental pressure signal are used to evaluate . This duration should also be changed according to the aims of the calibration and of the model use, since the design of a simple system could require the simulation of the first characteristic time, while the long-term analysis could require longer durations, with respect to that used here.

CONCLUSIONS

The implementation of the rheology of a polymeric material pipe in the transient governing equations requires the estimation of the viscoelastic parameters. The most common way to define such parameters is by means of a calibration procedure on experimental data. The main contribution of this paper is to address the problem of overfitting of viscoelastic models for transient analysis.

The implementation of the viscoelastic component in a frequency domain model reduces the computational time and speeds up the calibration. Furthermore, the presence of the viscoelastic parameters in the governing equations is limited to a complex function, . The analysis of this function for the considered calibrations shows that it could be used to enhance the differences from one model to another.

The standard procedure used, based on the minimization of the sum of the squares of the residuals by means of the combined use of a genetic and the Nelder–Mead algorithm, allows it to be shown that the Akaike criterion can be used to define the optimal number of parameters.

Since all the used algorithms and methodologies are easily available, this criterion, or other available criteria, can be easily used to avoid the overfitting. Additional investigations are needed to confirm some of the results shown and to asses the proper choice regarding different or more sophisticated viscoelastic models, optimization functions, calibration algorithms, and criteria for the comparison of different models.

REFERENCES

REFERENCES
Beaumont
,
M. A.
,
Zhang
,
W.
&
Balding
,
D. J.
2002
Approximate Bayesian computation in population genetics
.
Genetics
162
(
4
),
2025
2035
.
Burnham
,
K. P.
&
Anderson
,
D.
2003
Model Selection and Multimodel Inference: A Practical Information-Theoretic Approach
.
Springer
,
New York
.
Capponi
,
C.
,
Ferrante
,
M.
,
Zecchin
,
A. C.
&
Gong
,
J.
2017a
Experimental validation of the admittance matrix method on a y-system
.
Journal of Hydraulic Research
.
doi:10.1080/00221686.2017.1372818
.
Capponi
,
C.
,
Ferrante
,
M.
,
Zecchin
,
A. C.
&
Gong
,
J.
2017b
Leak detection in a branched system by inverse transient analysis with the admittance matrix method
.
Water Resources Management
10.1007/s11269-017-1730-6
.
Capponi
,
C.
,
Zecchin
,
A. C.
,
Ferrante
,
M.
&
Gong
,
J.
2017c
Numerical study on accuracy of frequency-domain modeling of transients
.
Journal of Hydraulic Research
.
doi:10.1080/00221686.2017.1335654
.
Chaudhry
,
M. H.
2014
Applied Hydraulic Transients
,
3rd edn
.
Springer
,
New York
.
Covas
,
D. I. C.
,
Stoianov
,
I.
,
Ramos
,
H.
,
Graham
,
N.
,
Maksimovic
,
C.
&
Butler
,
D.
2004
Water hammer in pressurized polyethylene pipes: conceptual model and experimental analysis
.
Urban Water Journal
1
(
2
),
177
197
.
Covas
,
D. I. C.
,
Stoianov
,
I.
,
Mano
,
J.
,
Ramos
,
H.
,
Graham
,
N.
&
Maksimovic
,
C.
2005
The dynamic effect of pipe-wall viscoelasticity in hydraulic transients. Part II – model development, calibration and verification
.
Journal of Hydraulic Research
43
(
1
),
56
70
.
Duan
,
H.-F.
,
Ghidaoui
,
M.
&
Lee
,
P. J.
2010a
Unsteady friction and visco-elasticity in pipe fluid transients
.
Journal of Hydraulic Research
48
(
3
),
354
362
.
Duan
,
H.-F.
,
Lee
,
P.
,
Ghidaoui
,
M. S.
&
Tung
,
Y.-K.
2010b
Essential system response information for transient-based leak detection methods
.
Journal of Hydraulic Research
48
(
5
),
650
657
.
Duan
,
H.-F.
,
Lee
,
P. J.
,
Ghidaoui
,
M. S.
&
Tung
,
Y.-K.
2012
System response function-based leak detection in viscoelastic pipelines
.
Journal of Hydraulic Engineering
138
(
2
),
143
153
.
Eaton
,
J. W.
,
Bateman
,
D.
,
Hauberg
,
S.
&
Wehbring
,
R.
2015
GNU Octave Version 4.0.0 Manual: a High-Level Interactive Language for Numerical Computations
.
Create Space Independent Publishing Platform
, .
Evangelista
,
S.
,
Leopardi
,
A.
,
Pignatelli
,
R.
&
de Marinis
,
G.
2015
Hydraulic transients in viscoelastic branched pipelines
.
Journal of Hydraulic Engineering
141
(
8
),
04015016-9
.
Ferrante
,
M.
&
Brunone
,
B.
2003
Pipe system diagnosis and leak detection by unsteady-state tests. 1. Harmonic analysis
.
Advances in Water Resources
26
(
1
),
95
105
.
Ferrante
,
M.
&
Capponi
,
C.
2017a
Experimental characterization of PVC-O pipes for transient modeling
.
Journal of Water Supply: Research and Technology – AQUA
66
(
8
),
606
620
.
doi:10.2166/aqua.2017.060
.
Ferrante
,
M.
&
Capponi
,
C.
2017b
Viscoelastic models for the simulation of transients in polymeric pipes
.
Journal of Hydraulic Research
.
doi:10.1080/00221686.2017.1354935
.
Ferrante
,
M.
,
Capponi
,
C.
,
Brunone
,
B.
&
Meniconi
,
S.
2015
Hydraulic characterization of PVC-O pipes by means of transient tests
.
Procedia Engineering
119
,
263
269
.
Ghilardi
,
P.
&
Paoletti
,
A.
1986
Additional viscoelastic pipes as pressure surges suppressors
. In:
5th International Conference on Pressure Surges, Hannover, Germany
(
Papworth
,
M.
, ed.).
BHRA, Cranfield
,
UK
, pp.
113
121
.
Ghilardi
,
P.
&
Paoletti
,
A.
1987
Viscoelastic parameters for the simulation of hydraulic transients in polymeric pipes
.
Excerpta of the Italian Contributions to the Field of Hydraulic Engineering
2
,
51
62
.
Gong
,
J.
,
Simpson
,
A. R.
,
Lambert
,
M. F.
,
Zecchin
,
A. C.
,
Kim
,
Y.-i.
&
Tijsseling
,
A. S.
2013
Detection of distributed deterioration in single pipes using transient reflections
.
Journal of Pipeline Systems Engineering and Practice
4
(
1
),
32
40
.
Jung
,
B. S.
&
Karney
,
B. W.
2008
Systematic exploration of pipeline network calibration using transients
.
Journal of Hydraulic Research
46
(
Suppl. 1
),
129
137
.
Kapelan
,
Z. S.
,
Savic
,
D. A.
&
Walters
,
G. A.
2003
A hybrid inverse transient model for leakage detection and roughness calibration in pipe networks
.
Journal of Hydraulic Research
41
(
5
),
481
492
.
Keramat
,
A.
,
Tijsseling
,
A. S.
,
Hou
,
Q.
&
Ahmadi
,
A.
2012
Fluid–structure interaction with pipe-wall viscoelasticity during water hammer
.
Journal of Fluids and Structures
28
(
C
),
434
455
.
Laucelli
,
D.
,
Romano
,
M.
,
Savić
,
D.
&
Giustolisi
,
O.
2016
Detecting anomalies in water distribution networks using EPR modelling paradigm
.
Journal of Hydroinformatics
18
(
3
),
409
427
.
Lee
,
P. J.
&
Vítkovský
,
J. P.
2010
Quantifying linearization error when modeling fluid pipeline transients using the frequency response method
.
Journal of Hydraulic Engineering
136
(
10
),
831
836
.
Lee
,
I.-Y.
,
Choi
,
S.-R.
&
Kang
,
M.-G.
2009
A method for measuring frequency series wave speed in viscoelastic pipes
. In:
IEEE Sensors Conference
,
Christchurch
,
New Zealand
, pp.
497
501
.
Lee
,
P. J.
,
Duan
,
H.-F.
,
Ghidaoui
,
M.
&
Karney
,
B.
2013
Frequency domain analysis of pipe fluid transient behaviour
.
Journal of Hydraulic Research
51
(
6
),
609
622
.
Massari
,
C.
,
Yeh
,
T.-C. J.
,
Ferrante
,
M.
,
Brunone
,
B.
&
Meniconi
,
S.
2014
Detection and sizing of extended partial blockages in pipelines by means of a stochastic successive linear estimator
.
Journal of Hydroinformatics
16
(
2
),
248
258
.
Meniconi
,
S.
,
Brunone
,
B.
,
Ferrante
,
M.
&
Massari
,
C.
2011
Transient tests for locating and sizing illegal branches in pipe systems
.
Journal of Hydroinformatics
13
(
3
),
334
345
.
Nash
,
J. E.
&
Sutcliffe
,
J. V.
1970
River flow forecasting through conceptual models. Part I – A discussion of principles
.
Journal of Hydrology
11
(
2
),
109
128
.
Pezzinga
,
G.
2000
Evaluation of unsteady flow resistances by quasi-2D or 1D models
.
Journal of Hydraulic Engineering
126
(
10
),
778
785
.
Pezzinga
,
G.
2002
Unsteady flow in hydraulic networks with polymeric additional pipe
.
Journal of Hydraulic Engineering
128
(
2
),
238
244
.
Pezzinga
,
G.
2014
Evaluation of time evolution of mechanical parameters of polymeric pipes by unsteady flow runs
.
Journal of Hydraulic Engineering
140
(
12
),
04014057-1–10
.
Pezzinga
,
G.
&
Scandura
,
P.
1995
Unsteady flow in installations with polymeric additional pipe
.
Journal of Hydraulic Engineering
121
(
11
),
802
811
.
Pezzinga
,
G.
,
Brunone
,
B.
,
Cannizzaro
,
D.
,
Ferrante
,
M.
,
Meniconi
,
S.
&
Berni
,
A.
2014
Two-dimensional features of viscoelastic models of pipe transients
.
Journal of Hydraulic Engineering
140
(
8
),
04014036-1
04014036-9
.
Sadegh
,
M.
&
Vrugt
,
J. A.
2014
Approximate Bbayesian computation using Markov chain Monte Carlo simulation: DREAM (ABC)
.
Water Resources Research
50
(
8
),
6767
6787
.
Savic
,
D. A.
,
Kapelan
,
Z. S.
&
Jonkergouw
,
P. M. R.
2009
Quo vadis water distribution model calibration?
Urban Water Journal
6
(
1
),
3
22
.
Soares
,
A. K.
,
Covas
,
D. I. C.
&
Reis
,
L. F. R.
2011
Leak detection by inverse transient analysis in an experimental PVC pipe system
.
Journal of Hydroinformatics
13
(
2
),
153
166
.
Sun
,
J.
,
Wang
,
R.
&
Duan
,
H.-F.
2016
Multiple-fault detection in water pipelines using transient-based time-frequency analysis
.
Journal of Hydroinformatics
18
(
6
),
975
989
.
Vítkovský
,
J. P.
,
Lee
,
P.
,
Zecchin
,
A. C.
,
Simpson
,
A. R.
&
Lambert
,
M. F.
2011
Head- and flow-based formulations for frequency domain analysis of fluid transients in arbitrary pipe networks
.
Journal of Hydraulic Engineering
137
(
5
),
556
568
.
Wylie
,
E. B.
&
Streeter
,
V.
1993
Fluid Transients in Systems
.
Prentice-Hall
,
Eaglewood Cliffs, NJ
.