Minimization of the Settling Time of Variable Area Flowmeters

In the article a differential equation describing transient behavior of variable area (VA) meters has been developed and validated experimentally for air as a measured fluid and for two float shapes—plumb bob and sphere. A modified version of simplex algorithm adapted for nonlinear constraint optimization problems was applied to minimize the settling time of VA meters in two cases. In the first case both the float and tube geometry were altered. In the second case only the float geometry was modified. The second case has been validated experimentally. The theory and experiment is in reasonable agreement (under 5% of full scale), which is satisfactory for the purposes of optimization of VA flowmeters dynamic performance. Analytical model of VA flowmeter has been proven to be a proper tool for optimization. Settling times obtained during the optimization process were several times shorter than these of commercially manufactured instruments. Overshoot has not exceeded the assumed value of 3%.


Introduction
Over the recent decades various versions of variable area (VA) flowmeters have grown to maturity. E.g., classic "see through" rotameters, orifice VA meters or VA flow transducers. Both the calibration procedures [1,2] and corrections for density and viscosity [3,4] have been established for a long time. Many float types have been thoroughly investigated [5].
The simplest indicating device is a transparent tube. More sophisticated flow transducers use a magnetically inert metal tube. In this case the float is provided with a permanent magnet and float position is sensed by an external magnet. Float displacement can then be converted to any standardized signal. Simple solution consists of a linear variable differential transformer [6]. Another possibility is the application of Hall sensors [7,8], magnetoresistive sensors matrix parallel to the tube [9] or capacitance transmitters [10]. In this case transducer does not exert forces on the float and therefore does not influence the VA meter dynamics. This simple and economical flowmeter continues to be one of the most widely used instruments for low and medium flow rates in numerous industry branches [11,12]. Industry 4.0, which is the current trend of automation and data exchange in manufacturing technologies including cyber-physical systems and Internet of Things, favors measurement sensors which have good dynamic performance [13], e.g., Coriolis flowmeters. However, price of Coriolis flowmeters is one or two orders of magnitude higher than price of simple VA flow transmitters. Study on dynamic performance of VA meters would allow optimization and development of simple and rugged flow transmitters that can meet contemporary needs of process automation.
Up till now research on dynamics of variable area flowmeters has concerned their behavior in presence of pulsations [14][15][16]. This aspect is important in case of such pulsation sources as reciprocating pumps or other sources classified in Ref. [17]. However, for control purposes such parameters as settling time, nature of step response (periodic or aperiodic) and overshoot are vital. These subjects have not yet been studied in the context of VA meters. The paper deals with these aspects and presents optimization procedures. An analytical model describing VA meter dynamics is developed and validated. The methodology bases on registering float position by a camera. This is the only method that does not exert forces on the float. The proposed equation has been proven to work for air and for two float shapes-plumb bob and sphere. Further research is required to validate the developed equation for other fluids and other float shapes. Computational fluid dynamics may be a viable tool for this purpose. In that case simulation validation criteria described in Ref. [18] may be applied. After validation of the transient equation the article present results of optimization for two cases. In the first case both the tube and float geometry were altered, in the second case only the float geometry was modified. The second case was validated experimentally and the presented methodology has been proven to be of scientific and practical significance.

VA Flowmeter Characteristics in Steady State
Variable area flowmeter consists of a tapered tube and a float. The float can take various forms. Examples are traditional plumb bob (as in Figure 1a), sphere or flat disc. Let us consider two perpendicular sections of the flowmeter marked in the Figure 1b. For stationary, incompressible flow one can use Bernoulli's law to write: where ρ is the fluid density and p 1 , U 1 , p 2 , U 2 are respectively pressure and mean velocities in Section 2.1. and Section 2.2.
Sensors 2019, 19 FOR PEER REVIEW 2 Up till now research on dynamics of variable area flowmeters has concerned their behavior in presence of pulsations [14][15][16]. This aspect is important in case of such pulsation sources as reciprocating pumps or other sources classified in Ref. [17]. However, for control purposes such parameters as settling time, nature of step response (periodic or aperiodic) and overshoot are vital. These subjects have not yet been studied in the context of VA meters. The paper deals with these aspects and presents optimization procedures. An analytical model describing VA meter dynamics is developed and validated. The methodology bases on registering float position by a camera. This is the only method that does not exert forces on the float. The proposed equation has been proven to work for air and for two float shapes-plumb bob and sphere. Further research is required to validate the developed equation for other fluids and other float shapes. Computational fluid dynamics may be a viable tool for this purpose. In that case simulation validation criteria described in Ref. [18] may be applied. After validation of the transient equation the article present results of optimization for two cases. In the first case both the tube and float geometry were altered, in the second case only the float geometry was modified. The second case was validated experimentally and the presented methodology has been proven to be of scientific and practical significance.

VA flowmeter Characteristics in Steady State
Variable area flowmeter consists of a tapered tube and a float. The float can take various forms. Examples are traditional plumb bob (as in Figure 1a), sphere or flat disc. Let us consider two perpendicular sections of the flowmeter marked in the Figure 1b. For stationary, incompressible flow one can use Bernoulli's law to write: where ρ is the fluid density and p1, U1, p2, U2 are respectively pressure and mean velocities in Sections 1.1. and 2.2.  The area of the narrowest part of the stream A 2 equals: and Cc is the contraction coefficient. As the velocity U 1 is much smaller than U 2 , its square can be neglected. Therefore: this pressure difference multiplied by the area of the float A f = πd f 2 /4 results in active force directed upwards: The restoring force F R is equal to the apparent weight of the float of volume V f and density ρ f immersed in the fluid of density ρ: where g-acceleration due to gravity. In the equilibrium state active force equals restoring force, therefore one can write: There are various methods to obtain Cc coefficient. Procedures of VA meter calculations based on d:d f diameter ratio and Ruppel number are shown in Ref. [2]. Another option is to determine Cc during calibration, as proposed by Ref. [1]. In this way a flow coefficient C(Re) is obtained, instead of the contraction coefficient Cc. Although d:d f is an important factor, it is taken here into account indirectly, through hydraulic diameter in Re. The flow coefficient C(Re) includes not only the effect of contraction, but also viscous friction, dynamic drag and other effects. In this case Re is defined as: where (d − d f ) is hydraulic diameter, U = Q/A(h)-mean velocity in the annular slot between float and tube and ν-kinematic viscosity of the fluid.

VA Behavior During Transient
During transient the active force is proportional to the velocity relative to the float. Therefore, float velocity dh/dt must be subtracted from velocity U 2 given by Equation (3). What is more, moving float acts as a piston pump. It results in an additional volumetric flow rate Q equal to: The resultant velocity equals: taking the above into the account, mean velocity relative to the float in the Section 2.2. is: the flow rate in the Section 2.2. equals: (12) and the active force during transient is then: taking into account dynamics, the restoring force equals: as Q can be in general a function of time Q(t), the equation of motion of the float will be therefore: It is easy to validate that in steady state (h = const.) this equation takes the form of static Equation (7). The relationship between the area A and the position of the float in the tube h (which is in fact the measure of the flow rate Q) is as follows: during transient Reynolds number equals: where: which allows to obtain the final equation: The differential Equation (19) is nonlinear and it can be only solved numerically. Runge-Kutta algorithm has been used for this purpose. There was a question of approximating function C(Re), which has been obtained experimentally. After numerous attempts with polynomial, exponential and logarithmic curves the best fit was obtained using approximating function in form of: where D and E coefficients were determined with the use of least squares method. As all the experiments were conducted with the use of atmospheric air as working medium, buoyancy forces could have been neglected. Therefore:

Validation of the Transient Equation
Seven VA flowmeters were investigated (Table 1). d low and d high are tube diameters at the heights corresponding to the lower and upper flow rate limits. ∆h is the distance between these extreme float positions (see Figure 1a). The tubes for VA flowmeter no. 6 and 7 are the same as for no. 2 and 4, only the float was replaced. The experiments were conducted with a methodology fully developed by the authors. The step response was filmed with an analog 8 mm camera with frame rate of 50 frames per second ( Figure 2a). The position of the float was later measured on the screen of a film previewer.
Step change in the flow rate was generated by a fast opening of a spring operated valve. The step was not ideal, therefore it has been registered by means of hot wire anemometer [19] and later approximated by a function in form of: where QSS is the flow rate in steady state. This form of excitation signal was used during the solution of the Equation (21). Exemplary signal registered by thermoanemometer and its approximation is The step was not ideal, therefore it has been registered by means of hot wire anemometer [19] and later approximated by a function in form of: where Q SS is the flow rate in steady state. This form of excitation signal was used during the solution of the Equation (21). Exemplary signal registered by thermoanemometer and its approximation is shown in the Figure 3. Test stand for settling time measurement is presented in the Figure 2b. The movement of the float at the final phase was slow. For that reason the overshoot could have been easily assessed by visual observation.
(a) (b) Figure 2. (a) Test stand for registration of the step response of VA flowmeter. 1-fast opening valve, 2-flow control valve, 3-constant temperature anemometer with data acquisition system, 4-VA flowmeter, 5-camera; (b) test stand for settling time measurement. 1-fast opening valve, 2-flow control valve, 3-VA flowmeter, 4-photoelectric sensors. Photoelectric sensors triggered a timer, enabling the measurement of the settling time.
The step was not ideal, therefore it has been registered by means of hot wire anemometer [19] and later approximated by a function in form of: where QSS is the flow rate in steady state. This form of excitation signal was used during the solution of the Equation (21). Exemplary signal registered by thermoanemometer and its approximation is shown in the Figure 3. Test stand for settling time measurement is presented in the Figure 2b. The movement of the float at the final phase was slow. For that reason the overshoot could have been easily assessed by visual observation. Step change registered by hot wire anemometer (solid line) and its approximation (dashed line).

Optimization
For a given flow rate range and measured fluid there are three parameters that define remaining parameters of a given VA flowmeter unambiguously. Examples are dlow, dhigh, Δh or mf, dlow, Δh. A set Step change registered by hot wire anemometer (solid line) and its approximation (dashed line).

Optimization
For a given flow rate range and measured fluid there are three parameters that define remaining parameters of a given VA flowmeter unambiguously. Examples are d low , d high , ∆h or m f , d low , ∆h. A set of parameters d f , (d low − d f ), ∆h was chosen for the reason of simplicity of formulae to calculate the remaining parameters. A modified simplex algorithm described in the Appendix A was used for optimization. Considering the parameters of VA flowmeters manufactured commercially, constraints presented in the Table 2 were assumed. The first three constraints are boundary constraints. As the density of float ρ f does not occur explicitly in the Equation (21), the 4 th constraint is referred to as inequality constraint. Volume of the float of a plumb bob type is expressed as V f = 0.679d 3 f . To determine values of coefficients D and E of the function approximating C(Re), all the experimental data concerning VA meters from Table 1 was collected. Function in a form of: was determined with the use of least squares method and was used during optimization. From the point of view of automation systems, minimum settling time and an aperiodic response (possibly with one small overshoot) is advantageous. The objective function was therefore defined as in Figure 4. The first three constraints are boundary constraints. As the density of float ρf does not occur explicitly in the Equation (21), the 4 th constraint is referred to as inequality constraint. Volume of the float of a plumb bob type is expressed as Vf = 0.679d 3 f.
To determine values of coefficients D and E of the function approximating C(Re), all the experimental data concerning VA meters from Table 1 was collected. Function in a form of: was determined with the use of least squares method and was used during optimization. From the point of view of automation systems, minimum settling time and an aperiodic response (possibly with one small overshoot) is advantageous. The objective function was therefore defined as in Figure  4. Optimized VA meters should have had a minimal settling time t0.97. The overshoot should have not exceeded 3% of the float position in a steady state hSS. For the each set of chosen VA parameters generated randomly by optimization algorithm described in the appendix, the values of an objective function were obtained by integration of the Equation (21). After every step of integration (ca. 0.01 s) an analysis of the step response was performed. If the float height exceeded the value of 1.03hSS (point a) the integration was interrupted as not fulfilling the assumptions. The results were then rejected and a new integration process was started, based on a new set of parameters generated by the optimization algorithm. If the response did not exceed the 1.03hSS value and in subsequent integration step the calculated h value started to decrease (point b), it was no longer possible that the float went beyond the ±0.03hSS limits. The value t0.97 calculated as an abscissa of the intersection of the line connecting points (ti, hi) and (ti+1, hi+1) and the horizontal line h = 0.97hSS was stored for further processing by the algorithm. To save time, if after For the each set of chosen VA parameters generated randomly by optimization algorithm described in the appendix, the values of an objective function were obtained by integration of the Equation (21). After every step of integration (ca. 0.01 s) an analysis of the step response was performed. If the float height exceeded the value of 1.03h SS (point a) the integration was interrupted as not fulfilling the assumptions. The results were then rejected and a new integration process was started, based on a new set of parameters generated by the optimization algorithm. If the response did not exceed the 1.03h SS value and in subsequent integration step the calculated h value started to decrease (point b), it was no longer possible that the float went beyond the ±0.03h SS limits. The value t 0.97 calculated as an abscissa of the intersection of the line connecting points (t i , h i ) and (t i+1 , h i+1 ) and the horizontal line h = 0.97h SS was stored for further processing by the algorithm. To save time, if after some arbitrarily chosen time the response did not reach 0.97h SS , the integration was interrupted (point c) and results were rejected. The optimization procedure stopped when five optimization results (objective function values) remained within a specified error band.

Validation of the Transient Equation
The comparison of experiments with the results obtained from the Equation (21) are presented in Figures 5-10. One can distinguish between three typical step responses. For VA no. 1 an aperiodic response was obtained (Figures 5 and 6). It resembles a step response of a first order linear system with damping coefficient above 1. The settling time rises with the decreasing value of step change, which is in contrast to a typical linear system. The comparison of experiments with the results obtained from the Equation (21) are presented in Figures 5-10. One can distinguish between three typical step responses. For VA no. 1 an aperiodic response was obtained (Figures 5 and 6). It resembles a step response of a first order linear system with damping coefficient above 1. The settling time rises with the decreasing value of step change, which is in contrast to a typical linear system.  For VA no. 5 results resemble a step response of a second order linear system with a small value of damping coefficient ( Figure 7). As opposed to the linear system, the oscillation frequency is higher for a lower step value. Moreover, oscillation fragments below steady-state value are shorter than the upper parts. The comparison of experiments with the results obtained from the Equation (21) are presented in Figures 5-10. One can distinguish between three typical step responses. For VA no. 1 an aperiodic response was obtained (Figures 5 and 6). It resembles a step response of a first order linear system with damping coefficient above 1. The settling time rises with the decreasing value of step change, which is in contrast to a typical linear system.  For VA no. 5 results resemble a step response of a second order linear system with a small value of damping coefficient ( Figure 7). As opposed to the linear system, the oscillation frequency is higher for a lower step value. Moreover, oscillation fragments below steady-state value are shorter than the upper parts. The response of VA no. 3 resembles a step response of a linear system with damping coefficient equal to ca. 0.5, practically with one overshoot (Figure 8). In contrast to the linear system, the overshoot decreases with a decreasing step value.  The response of VA no. 3 resembles a step response of a linear system with damping coefficient equal to ca. 0.5, practically with one overshoot (Figure 8). In contrast to the linear system, the overshoot decreases with a decreasing step value. The response of VA no. 7 of a sphere shape (Figures 9 and 10) is similar to the response of VA no. 5, but with a smaller value of damping coefficient.  Uncertainty analysis is shown in the Appendix B. The sum of the experimental and theoretical uncertainties is in the worst case equal to U95(h) = 18.2 mm whereas the differences between experiments and simulations were up to 11 mm. It was assumed that the mathematical model was in a satisfactory agreement with the experiments.  Uncertainty analysis is shown in the Appendix B. The sum of the experimental and theoretical uncertainties is in the worst case equal to U95(h) = 18.2 mm whereas the differences between experiments and simulations were up to 11 mm. It was assumed that the mathematical model was in a satisfactory agreement with the experiments. For VA no. 5 results resemble a step response of a second order linear system with a small value of damping coefficient ( Figure 7). As opposed to the linear system, the oscillation frequency is higher for a lower step value. Moreover, oscillation fragments below steady-state value are shorter than the upper parts.

Optimization of VA Flowmeters
The response of VA no. 3 resembles a step response of a linear system with damping coefficient equal to ca. 0.5, practically with one overshoot (Figure 8). In contrast to the linear system, the overshoot decreases with a decreasing step value.
The response of VA no. 7 of a sphere shape (Figures 9 and 10) is similar to the response of VA no. 5, but with a smaller value of damping coefficient.
Uncertainty analysis is shown in the Appendix B. The sum of the experimental and theoretical uncertainties is in the worst case equal to U 95 (h) = 18.2 mm whereas the differences between experiments and simulations were up to 11 mm. It was assumed that the mathematical model was in a satisfactory agreement with the experiments.

Optimization of VA Flowmeters
Three VA meters of a plumb bob type were optimized for nominal ranges listed in the Table 3. In studied cases the pipe length ∆h was equal to the lower bound l i = 150 mm and the overshoot was equal to 3% of h SS value. In comparison to the commercial VA meters, resultant float density was relatively low and the clearance between a float and a tube was relatively high. The settling time (0.32-0.35) s was shorter than for commercial VA meters (see  where it equals (0.9-1.5) s). Unfortunately, unit production of glass tubes for experiments is costly. Therefore the optimization was repeated for the existing tubes specified in the Table 1, and was limited to modifying float parameters. In this case optimization algorithm generated the lower range limit Q min and the clearance (d low − d f ). The remaining parameters such as upper range limit, float mass and density were then easily calculated. Additional restriction was added concerning the rangeability, namely Q max /Q min ≥ 8.
Optimization results are presented in the Table 4. Columns 4 and 6 contain values of the clearance (d low − d f ) and Q min generated by the optimization algorithm. The calculated values of ρ f and Q max are presented in columns 5 and 7. Theoretical values of the settling time (column 8) are within the range of (0.341-0.372) s, which is also an improvement in comparison to the commercial meters. The experimental results are presented in columns 9-12. There is a good agreement between results of optimization and experiments. Measured settling times are equal or slightly shorter than theoretical values. The measured overshoot has never exceeded the assumed value of 3%.

Conclusions
Analytical model of VA flowmeters dynamics has been developed and validated experimentally. It has been proven to be a proper tool for optimization. This will make VA meters more appropriate for process automation, which is a prerequisite for the implementation of Industry 4.0 solutions. Settling times obtained during the optimization process were several times shorter than these of commercially manufactured instruments: (0.341-0.372) s instead of (0.9-1.5) s. This proves that the presented methodology is of scientific and practical significance. What is more, it has been proven that the developed equation describing VA meter dynamics is appropriate for at least two float shapes, i.e., plumb bob and sphere.
The sum of the experimental and theoretical uncertainties is in the worst case 18.2 mm whereas the differences between experiments and simulations were up to 11 mm. The discrepancies can also be attributed to the numerical modeling of the system (e.g., inaccurate physical models, discretization and round-off errors) and the fact that the uncertainty analysis has been developed for the steady state case, not for transient. What is more, the approximating function C(Re) that takes float shape into account has been obtained as the best fit of experimental data. Comparing the Figure 9 with Figure 5, Figure 7and Figure 8, it can be seen that the differences between modelling and experiments are smaller for a simple sphere float shape. However, the theory and experiment are in reasonable agreement (under 5% of full scale), satisfactory for the purposes of optimization of VA flowmeters dynamic performance. In generally, the discrepancies between the experiments and calculations are smaller for larger rotameters.
Further research should be conducted for various float shapes and various fluids as a measured medium. This includes both classical rotameters and variable orifice meters. The influence of forces from transmitter acting on a float should also be studied. It would also be beneficial to employ computational fluid dynamics for further optimization.

Appendix A. Description of the Optimization Algorithm
The algorithm used in the article is a modified version of a simplex algorithm adapted for nonlinear constraint optimization problems. A set of parameters (ω 1 , ω 2 , . . . , ω n ) that assures minimum of the objective function defined in Figure 4 is sought. Calculation of the objective function procedure is an external procedure from the point of view of the optimization algorithm.
At first five points (vertices) are randomly generated with the use of a formula: where r i is a random number from the range 0, i) l i , u i are the lower and upper bounds for the optimization variable. If the generated set does not fulfill all of the constraints, a new set of points is generated. The applied procedure operates on k ≥ n + 5 points. Subsequent points are also generated with the use of a random numbers generator. If the generated point does not fulfill the functional constraints, it is shifted half a distance towards the center of gravity (CG) of the existing set of points, e.g., from point 9 to 9 . Such a procedure is repeated till the point is located in the acceptable area.
The objective function is calculated in each node of the polygon formed by all the randomly generated points. After generating k points the point with the maximum (the worst) value of the objective function (e.g., point 5) is replaced by an another one. This point is placed on a straight line crossing this worst point and the CG, on the opposite side of CG. The distance between the new point (e.g., 10 in the Figure A1) and CG is equal to the distance between the worst point and CG multiplied by a reflection coefficient α > 1. optimization variable. If the generated set does not fulfill all of the constraints, a new set of points is generated.
The applied procedure operates on 5 k n ≥ + points. Subsequent points are also generated with the use of a random numbers generator. If the generated point does not fulfill the functional constraints, it is shifted half a distance towards the center of gravity (CG) of the existing set of points, e.g., from point 9 to 9′. Such a procedure is repeated till the point is located in the acceptable area. The objective function is calculated in each node of the polygon formed by all the randomly generated points. After generating k points the point with the maximum (the worst) value of the objective function (e.g., point 5) is replaced by an another one. This point is placed on a straight line crossing this worst point and the CG, on the opposite side of CG. The distance between the new point (e.g., 10 in the Figure A1) and CG is equal to the distance between the worst point and CG multiplied by a reflection coefficient α > 1.
If this point does not fulfill the hard constraints, it is shifted towards the acceptable area, at some distance ε from the border. If the point does not fulfill the functional constraints, it is shifted towards the CG by half the distance between this point and the CG. This procedure is repeated till the point reaches an acceptable area. Figure A1. Scheme of the optimization procedure.
If this point does not fulfill the hard constraints, it is shifted towards the acceptable area, at some distance ε from the border. If the point does not fulfill the functional constraints, it is shifted towards the CG by half the distance between this point and the CG. This procedure is repeated till the point reaches an acceptable area.
Introduction of the reflection coefficient α > 1 causes continuous growth of the area covered by simplex. It compensates the reduction of this area caused by the procedure of halving the distance between point lying outside the acceptable area and the GG. It makes it possible to find the minimum of the objective function even if the start points lie at a great distance from the optimal point. This way simplex does not collapse into a subspace.
The optimization procedure stops when five objective function values remain within a specified error band. After the tests of the algorithm it has turned out that k should lie in the range of (15)(16)(17)(18)(19)(20)(21)(22)(23)(24)(25). For k less than 10 it has been difficult to find the optimum and for k > 25 the calculation time has been unnecessarily prolonged.

Appendix B.1. Uncertainty of Calculations
An assessment of uncertainty of dynamic phenomena would be very difficult, if not impossible. Therefore, the analysis of uncertainty for steady state has been developed under assumption that the results will be comparable to transient.
According to Ref. [20,21], if the functional relationship between the input quantities X 1 , X 2 , . . . , X N and output quantity Y in a measurement process is: then an estimate of Y denoted by y is obtained from Equation (A2) using input estimates x 1 , x 2 , . . . , x N : For uncorrelated input quantities X i the combined uncertainty may be found by combining the uncertainty of each of the contributing factors: where u(x i ) are individual components of uncertainty and c i are sensitivity coefficients.
In VA meter analysis the output quantity is float elevation h. This applies both to the analytical and the experimental case. The basic functional relationship is a slightly transformed Equation (7): where: Sensitivity coefficients c i are defined as the ratio of change of the output quantity y with respect to the input quantity x i . If an explicit formula for the output quantity is given, values of sensitivity coefficients are obtained by partial differentiation. In the analyzed case an another method proposed in Ref. [21] is applied. Sensitivity coefficients are obtained numerically, by calculating the effect of a small change in the input variable on the output. Output value y is first calculated using x i , and then recalculated using (x i + ∆x i ), where ∆x i is a small increment in x i . The result can be written as (y + ∆y), where ∆y is the increment in y caused by ∆x i . Sensitivity coefficients are then calculated as: A functional relationship between the input quantities and the output quantity h can be written, taking into account Equations (A5) and (A6), as: f (h) = C(Re)π hd f tan γ + h 2 tan 2 γ 2gm Then the values of h can be calculated as the roots of the Equation (A8) with the use of standard numerical methods. Exemplary calculations of sensitivity coefficients for VA meter (3-30) dm 3 /h and flow rate step 10 dm 3 /h are presented in the Table A1. After the assessment of all the uncertainty components according to guidelines from Ref. [20], the uncertainty budget was established (Table A2). During the uncertainty assessment the following factors have been taken into account. Type B uncertainty evaluation was used based on manufacturers' declarations for d f (micrometer), d high and d low (stick micrometer), ∆h (caliper) and m (scales). For ρ, υ uncertainty of tabular data was used. For Q uncertainty of bell prover was used and for experimentally obtained C uncertainty was assessed using the uncertainty of measurement instruments and the goodness-of-fit of the approximating function (23). Finally, the combined uncertainty equals: (c i u(x i )) 2 = 9.01 mm (A9) and the expanded uncertainty is equal to: U 95 (h) = 2u c (h) = 18.022 mm (A10) Similar results were obtained for other flowmeters. The uncertainty was smaller for larger VA meters, e.g., for VA meter no. 5 for Q = 975 dm 3 /h U 95 (h) = 1.98 mm, for Q = 1500 dm 3 /h U 95 (h) = 2.77 mm and for VA meter no. 7 Q = 1381 dm 3 /h U 95 (h) = 4.46 mm.