Monte Carlo determination of the uncertainty of effective area and deformation coefficient for a piston cylinder unit


ACTA IMEKO | www.imeko.org December 2020 | Volume 9 | Number 5 | 338 

ACTA IMEKO 
ISSN: 2221-870X 
December 2020, Volume 9, Number 5, 338 - 342 

 
 

MONTE CARLO DETERMINATION OF THE UNCERTAINTY OF 

EFFECTIVE AREA AND DEFORMATION COEFFICIENT FOR A PISTON 

CYLINDER UNIT 
 

C. Wuethrich1, S. Souiyam2 

 
1 METAS, Bern, Switzerland, christian.wuethrich@metas.ch  

2 LPEE-LNM, Casablanca, Morocco, Souiyam@lpee.ma  

 

Abstract: 

This paper describes the approach made to 

describe the uncertainty on the area at zero pressure 

(A0) and on the deformation coefficient (𝜆) of 
piston-cylinders used for pressure definition. We 

perform the Monte-Carlo simulation for an ordinary 

least squares (OLS) as well as for a weighted least 

squares (WLS) and a generalized least squares 

(GLS). We also introduce an innovative technique 

to improve the results obtained by WLS in order to 

get close to the quality obtained using GLS. We 

discuss the different situations that guides to the best 

choice between OLS, WLS and GLS. 

Keywords: Monte-Carlo, uncertainty, pressure 

balance, piston-cylinder 

1. INTRODUCTION 

Traditionally, pressure laboratories have an 

uncertainty on the effective area determined by 

cross floatation based on the uncertainties on the 

influence factors and the sensitivity coefficients. [1-

3]. The uncertainty on the area at zero pressure (A0) 

and the deformation coefficient (𝜆) determined by 
least squares is not trivial to obtain by the technique 

of the sensitivity coefficients. The project Euramet 

1125 [4] described the results obtained by several 

national metrology institutes in Europe for a set of 

data, simulating measurements affected by known 

uncertainties. Alternative techniques have been 

recently proposed to improve the uncertainty [5-7]. 

We explore the uncertainty of the determination of 

effective area at zero pressure and deformation 

coefficient based on Monte-Carlo simulation of the 

calibration process and the determination of 𝜆 and 
A0 

2. GOAL OF THIS WORK 

METAS and LPEE-LNM performed recently a 

bilateral comparison of piston-cylinder units. We 

wanted not only to compare the effective area, at a 

given pressure, but also the values for A0 and 𝜆 

obtained by least squares calculation on the set of 

measurements. A Monte-Carlo simulation is an 

efficient and mathematically correct way to assess 

the uncertainty of a value obtained using a least 

squares calculation. We are confronted to 

inhomogeneous uncertainties on the effective area 

depending of the pressure at which the effective area 

is measured, leading to heteroscedasticity on the 

residuals. We also have a correlation in the error of 

the mass used on each pressure balance as the mass 

used for the realisation of the first pressure step 

remains on the pressure balance for the second 

pressure step and so on for each successive pressure 

step. The whole set of mass is also calibrated against 

the same reference and the different pieces of charge 

are correlated to some degree. 

 
Figure 1: Setup used for the determination of the piston-

cylinder effective area. The reference pressure balance, 

on the left, is working with gas and the pressure balance 

under calibration, on the right, is working with oil. 

3. DETERMINATION OF A0 AND LAMBDA 

In pressure metrology we determine the effective 

area of a piston cylinder unit (PCU) at different 

pressure points by cross floatation [2-3]. Based on 

these measurements we want to determine the 

effective area at zero pressure (A0) and the 

deformation coefficient (𝜆), in order to have a model 
of the deformation of the PCU under pressure. The 

area of the PCU at a given pressure is the given by: 

𝐴(𝑝) = 𝐴0(1 + 𝜆𝑝) (1) 

http://www.imeko.org/
mailto:christian.wuethrich@metas.ch
mailto:Souiyam@lpee.ma


ACTA IMEKO | www.imeko.org December 2020 | Volume 9 | Number 5 | 339 

Because we have only two free parameters to 

determine but have measurements performed at a 

larger number of pressure points it is convenient to 

use a least square approach to solve this equation. It 

is convenient to rewrite (1) the following way: 

𝐴(𝑝) = 𝐴0 + 𝑏 ∙ 𝑝 (2) 

Where b is in fact A0*𝜆. 

3.1. Least squares determination 
The model used for the relation between the set 

of measurement, the deformation coefficient and the 

effective area at zero pressure is as follow: 

(
1 𝑝1
⋮ ⋮

1 𝑝𝑛

) ∙ (
𝐴0
𝑏

) = (
𝐴(𝑝1)

⋮

𝐴(𝑝𝑛 )
) + (

𝜀1
⋮

𝜀𝑛
) . (3) 

Where the pi are the pressure steps at which we 

determined the effective area A(pi). The value A0 is 

the effective area at zero pressure and b is the 

increase of area per unit pressure due to deformation 

and these two values are the unknown of the system. 

The 𝜀i are the residual that we will minimise using 
the least square equation. 

We can rewrite the equation 2 using the 

following matrix definition: 

𝑋𝐶 = 𝑌 + 𝐸 (4) 

We introduce a weighting matrix that multiplies 

the terms of the equation on the left and on the right. 

𝑉 −1𝑋𝐶 = 𝑉 −1𝑌 (5) 

The matrix V is the weighing coefficient. V is an 

identity matrix in the case of ordinary least square 

(OLS). V is a diagonal matrix with coefficients 

proportional to the uncertainty of A(pi) or to the 

square of the uncertainty of A(pi) in the case of the 

weighted least squares according to respectively the 

uncertainty (WLS-U) or the uncertainty squared 

(WLS-U2). V is the matrix of variance and 

covariance of the A(pi) in the case of generalised 

least squares (GLS). For simplicity of writing we 

define the matrix W as V-1 

We then solve the equation by multiplying by X´, 

the transposate of X on the left side: 

𝑋′𝑉 −1𝑋𝐶 = 𝑋′𝑉 −1𝑌 (6) 

And finally we multiply by (X´V-1X)-1 on the left 

side and obtain the solution of the equation system: 

𝐶 = (𝑋′𝑉 −1𝑋)−1 ∗ 𝑋′𝑉 −1𝑌 (7) 

4. MONTE CARLO SIMULATION 

The Monte-Carlo simulation reproduces the 

errors obtained when performing the determination 

of the effective area of a piston-cylinder using the 

cross-floating technique. The effective area for a 

given pressure step is obtained by using the usual 

equation of effective area determination. In a first 

step we calculate the pressure generated by the 

reference pressure balance by the jth pressure step: 

𝑝𝑗 =
∑ 𝑚𝑖 𝑔(1 −

𝜌𝑎
𝜌𝑚𝑖⁄ )𝑖

𝐴𝑅0(1 − 𝜆𝑝) (1 + (𝛼𝑝 + 𝛼𝑐 )(𝑡 − 𝑡𝑟 ))
 . (8) 

And then we calculate the effective area of the 

piston under test. 

𝐴𝑝𝑗 =
∑ 𝑚𝑖𝑗 𝑔 (1 −

𝜌𝑎𝑗
𝜌𝑚𝑖𝑗⁄ )𝑖

𝑝𝑗 (1 + (𝛼𝑝 + 𝛼𝑐 )(𝑡 − 𝑡𝑟 ))
 . (9) 

In our model there is 9 pressure steps and we 

assume that the measurement is repeated 10 times 

in order to determine the effective area. A set of 

measurement comprises the 90 measurement 

resulting from the 9 pressure steps repeated ten 

times. 

We repeat the calculation to obtain 10'000 

measurement sets in order to perform the 

uncertainty calculation on A0 and 𝜆 of the piston 
under calibration. 

The errors we generate in the Monte-Carlo 

simulation are supposed to depict as close as 

possible the reality of the metrological process. All 

influence factors have an error contribution 

considered constant in the whole set of 

measurements that do not change from one pressure 

step to the next one. Some influence factors are 

considered to have an error contribution that does 

change at each pressure step. This is the case for 

influence factors that need a new measurement at 

each pressure step. 

 

 
Figure 2: Plot (blue stars) of the effective area obtained 

in the Monte-Carlo simulation. The red line is the average 

value of the effective area and the green curves are placed 

at plus and minus the standard deviation. 

 

4.1. Generation of the set of data 
As illustration of the numerical technique we 

simulate a piston of 2 bar/kg working with air used 

to calibrate a piston of 5 bar/kg working with oil. 

0e00 1e072e06 4e06 6e06 8e06 1.2e07 1.4e07 1.6e07

1.962e-05

1.9616e-05

1.9618e-05

1.9622e-05

1.9624e-05

1.9617e-05

1.9619e-05

1.9621e-05

1.9623e-05

Pressure (Pa)

A
re

a
 (

m
2
)

Area versus pressure used in the simulation

http://www.imeko.org/


ACTA IMEKO | www.imeko.org December 2020 | Volume 9 | Number 5 | 340 

We perform the calibration at 1 MPa, 2 MPa and 

then each 2 MPa up to 16 MPa. The values and 

respective uncertainties type A and B are given in 

table 1. We assume for all the uncertainties a normal 

distribution in our simulation. 

The simulation is performed using scilab 6.1.0 [9] 

on a small computer running a 64 bits operating 

system like Linux or Windows. The generation of 

the 10'000 sets of data, the resolution of the least 

squares equations for the different weighting 

schemes and some statistical calculations takes 

about 2 minutes. 

 

Table 1: Values and related uncertainties of type A and type B used in the simulation. 

Influence factor Abreviation Value Uncertainty 

type B 

Uncertainty 

type A 

Effective area reference piston A0,ref  4.905 × 10
-5 m2 9.64 × 10-10 m2  

Deformation coefficient reference 𝜆ref 7.0 × 10-13 Pa-1 1.0 × 10-13 Pa-1  

Effective area piston to calibrate A0,cal  1.962 × 10
-5 m2   

Def. coefficient piston to calibrate 𝜆cal  7.0 × 10
-13 Pa-1   

Thermal expansion coeff. reference αref 9 ppm/°C 1 ppm/°C  

Therm. Exp. coeff. Piston to calibrate αcal 9 ppm/°C 1 ppm/°C  

Mass used on the reference piston Mi,ref 5 kg 50 mg  

10 kg 50 mg  

Mass used on the piston to calibrate Mi,cal 2.0 kg 20 mg  

4.0 kg 20 mg  

Height of the oil column Δh 0 m 0.05 m 0.003 m 

Temperature of the reference piston Tref 22 °C 0.1 °C 0.1 °C 

Temperature of the piston to calibrate Tcal 22 °C 0.1 °C 0.1 °C 

Density of the mass used on the ref. ρMiref 7950 kg/m
3 71 kg/m3  

Dens. of the mass used on the sample ρMical 7950 kg/m
3 71 kg/m3  

Density of air ρair 1.2 kg/m
3 0 kg/m3 0.005 kg/m3 

 

5. DETERMINATION OF THE 
WEIGHTING MATRIX 

The weighting matrix V used in the least square 

calculation according to equation 5 can be defined 

of different manners [8]. The simplest way is to take 

a diagonal matrix with 1. This is the situation of 

ordinary least squares determination (OLS). 

𝑣𝑖𝑗 = 0: 𝑖 ≠ 𝑗 ,   𝑣𝑖𝑖 = 1  (10) 

Traditionally people working on the 

determination of the parameters of PCU use a 

diagonal weighting matrix in which the values of the 

elements on the diagonal are related to the 

uncertainty of the effective area obtained at a given 

pressure step. (WLS-U) 

𝑣𝑖𝑗 = 0: 𝑖 ≠ 𝑗 ,   𝑣𝑖𝑖 = 𝑢(𝐴(𝑝𝑖 )) (11) 

Another solution consists in a weighting matrix 

in which the diagonal elements are determined by 

the square of the uncertainty of the effective area at 

a given pressure. (WLS-U2) 

𝑣𝑖𝑗 = 0: 𝑖 ≠ 𝑗 ,   𝑣𝑖𝑖 = 𝑢
2(𝐴(𝑝𝑖 )) (12) 

Finally we can determine the matrix V using the 

variance and covariance of the effective area at the 

different pressure steps. This approach is 

susceptible to provide the best uncertainties but is 

difficult to implement because the determination of 

the covariance can be challenging. In our situation 

we have a large set of simulated data through the 

Monte-Carlo calculation and are able to calculate 

the covariance based on the set of data at our 

disposal. (GLS) 

𝑣𝑖𝑗 = 𝑐𝑜𝑣 (𝐴(𝑝𝑖 ), 𝐴(𝑝𝑗 )) : 𝑖 ≠ 𝑗 ,    

𝑣𝑖𝑖 = 𝑢
2(𝐴(𝑝𝑖 )) 

(13) 

Until now all the matrices we have considered 

are based on the uncertainties observed on the 

effective area of the PCU. In a recent publication [6], 

P. Otal introduced the notion of matrix coefficients 

based on the uncertainty of Ã(p) which is calculated 

assuming only the uncertainty contributions on the 

effective area that are pressure dependant and 

contribute to heteroscedasticity. In our model we 

take into account the contribution of the height of 

the oil column (Δh) and the deformation coefficient 

of the reference PCU (𝜆ref). The plot of the Ã(p) and 
their standard deviation is shown on Fig. 3. 

We are able, based on the uncertainty of Ã(p), to 

define the following matrices used in the least 

squares determination. 

http://www.imeko.org/


ACTA IMEKO | www.imeko.org December 2020 | Volume 9 | Number 5 | 341 

A weighted least squares based on the 

uncertainty of Ã(p) as explained in [6] (WLS-UÃ): 

𝑣𝑖𝑗 = 0: 𝑖 ≠ 𝑗 ,   𝑣𝑖𝑖 = 𝑢 (�̃�(𝑝𝑖 )) (14) 

A weighted least squares based on the square of 

the uncertainty of Ã(p) (WLS-U2Ã): 

𝑣𝑖𝑗 = 0: 𝑖 ≠ 𝑗 ,   𝑣𝑖𝑖 = 𝑢
2 (�̃�(𝑝𝑖 )) (15) 

Finally, we can also define a generalised least 

squares matrix based on Ã(p) (GLS-Ã): 

𝑣𝑖𝑗 = 𝑐𝑜𝑣 (�̃�(𝑝𝑖 ), �̃�(𝑝𝑗 )) : 𝑖 ≠ 𝑗 ,    

𝑣𝑖𝑖 = 𝑢
2 (�̃�(𝑝𝑖 )) 

(16) 

 

 
Figure 3: Plot (blue stars) of the effective area obtained 

in the Monte-Carlo simulation taking into account the 

uncertainties related to the oil column and the 

deformation coefficient. The red line is the average value 

of the effective area and the green curves are placed at 

plus and minus the standard deviation. 

It is interesting to note that the diagonal elements 

of the matrix are similar in the equation 12 and 13 

as well as in equation 15 and 16. 

The coefficients of the diagonal of the matrix W 

are summarised in Table 2 and show the relative 

importance given to the different measurements. 

Table 2: Elements wii used according to the different 

weighing schemes. A normalisation has been applied so 

that the larger element is equal to one. 

 WLS 

-U 

WLS 

-U2 

WLS 

-UÃ 

WLS 

-U2Ã 

w11 0.395 0.157 0.086 0.008 

w22 0.644 0.414 0.170 0.029 

w33 0.861 0.743 0.325 0.106 

w44 0.937 0.877 0.468 0.219 

w55 0.968 0.938 0.600 0.361 

w66 0.984 0.968 0.722 0.521 

w77 0.992 0.984 0.830 0.689 

w88 0.997 0.993 0.923 0.851 

w99 1.000 1.000 1.000 1.000 

 

5.1. Uncertainties on A0 and lambda 
We applied the different weighing schemes to 

our set of data. We managed to obtain an average 

value for A0 and 𝜆 similar to the nominal value for 
all the weighting matrices demonstrating that the 

technique is correct. The standard deviation 

obtained on A0 and 𝜆 greatly depends on the 
weighting scheme selected. The results summarised 

on table 3 show that the GLS techniques delivers the 

best uncertainties. Weighting least squares based on 

the square of the uncertainty are more performant 

than those based on the uncertainty. The most 

interesting result is that the uncertainty on A0 and 𝜆 
for the WLS based on u2(Ãp) is the closest from the 

results obtained with a GLS approach. It is also 

interesting to note that the GLS approach based on 

Ã(p) delivers results similar with the GLS based on 

A(p). 

Table 3: Values and respective uncertainties obtained for 

the area at zero pressure and the deformation coefficient 

for the PCU under calibration, for different least squares 

techniques. 

 A0 𝜆 
Nominal value 1.962 × 10-5 m2 0.70 × 10-12 Pa-1 

u(OLS) 36 ppm 2.20 × 10-12 Pa-1 
u(WLS-U) 31 ppm 1.70 × 10-12 Pa-1 

u(WLS-U2) 28 ppm 1.35 × 10-12 Pa-1 

u(GLS) 20 ppm 0.36 × 10-12 Pa-1 

u(WLS-UÃ) 28 ppm 1.21 × 10-12 Pa-1 

u(WLS-U2Ã) 24 ppm 0.74 × 10-12 Pa-1 

u(GLS-Ã) 20 ppm 0.35 × 10-12 Pa-1 

 

The plot of 𝜆 versus A0 for the different least 
square weighting matrices, depicted in Fig. 4, shows 

the uncertainties and correlation obtained with the 

different least squares approximations. 

In order to assess the independence of the 

uncertainties on A0 and 𝜆, we calculated the 
correlation between the two parameters for the 

different weighing matrices used in the least square 

calculation. The results, presented on table 4, show 

that the correlation is relatively high when OLS is 

used but can be decreased to almost negligible 

amount in the case of a GLS approach. The 

correlation for the WLS according to u2(Ãp) is the 

best value if we exclude the GLS. 

Table 4: Correlation between A0 and 𝜆 for the different 
weighing matrices tested in the simulation. 

Weighing matrix Correlation A0 - 𝜆 
OLS -0.83 
WLS-U -0.77 

WLS-U2 -0.71 

GLS -0.28 

WLS-UÃ -0.69 

WLS-U2Ã -0.55 

GLS-Ã -0.27 

 

0e00 1e072e06 4e06 6e06 8e06 1.2e07 1.4e07 1.6e07

1.962e-05

1.9618e-05

1.9622e-05

1.9617e-05

1.9619e-05

1.9621e-05

1.9623e-05

1.96175e-05

1.96185e-05

1.96195e-05

1.96205e-05

1.96215e-05

1.96225e-05

Pressure (Pa)

A
re

a
 (

m
2
)

Area versus pressure used in the simulation

http://www.imeko.org/


ACTA IMEKO | www.imeko.org December 2020 | Volume 9 | Number 5 | 342 

Figure 4: Plot of 𝜆 versus A0 for different least square calculation: From left to right, top:OLS, WLS-U , WLS-U2 and 
GLS, bottom: WLS_UÃ, WLS_U2Ã, GLS_Ã. 

 
 

6. SUMMARY 

We applied successfully a Monte-Carlo 

simulation to determine the uncertainty on A0 and 𝜆 
determined by cross floating of pressure balances. 

Our numerical simulation has been used to optimise 

the weighting matrix used in the determination of A0 

and 𝜆. The raw data from the Monte-Carlo 
simulation have been used to determine the 

covariance between the effective areas determined 

at two different pressure. 

We have shown that an improvement of the 

weighting matrix proposed by Otal and Al. is able 

to provide an uncertainty close to what is obtained 

using a generalised least squares approach. 
The authors want to thank the PTB for providing 

this collaboration opportunity between METAS and 

LPEE-LNM through the project "Promotion of 

quality-assurance capabilities and services in the 

Maghreb to strengthen international trade". We 

want to thank Dr. Kilian Marti for his help in 

correcting the manuscript. 

 

7. REFERENCES 

[1] Dadson R S, Lewis S L and Peggs G N, The Pressure 
Balance: Theory and Practice (London: HMSO) pp 

19–72, 1982. 

[2] RMAéro 802-21, Étalonnage et utilisation des 
balances manométriques—Balance à application de 

masses 

[3] Calibration of Pressure Balances, EURAMET cg-3, 
Version 1.0 (03/2011) EURAMET e.V., Technical 

Committee for Mass 

[4] I. Morgado et al, "Evaluation of effective float 
measurement with pressure balances", EURAMET 

project 1125, final report 

[5] Sutton C M, Fitzgerald M P and Giardini W, A 
method of analysis for cross-floats between pressure 

balances, Metrologia, vol. 42, pp. 212–215, 2005. 

[6] P. Otal, C. Yardin, Modelling methods for pressure 
balance calibration, Measurement Science and 

Technology, vol. 31, 034004, 2020. 

[7] V. Ramnath, Determination of pressure balance 
distortion coefficient and zero-pressure effective 

area uncertainties, Int. J. Metrol. Qual. Eng, vol. 2, 

pp. 101-119, 2011.  

[8] Strutz and Tilo, Data Fitting and Uncertainty (A 
practical introduction to weighted least squares and 

beyond)., 2010 

[9] https://www.scilab.org 
 

 

1.962e-051.9618e-05 1.9622e-05

0e00

-5e-12

5e-12

Area (m2)

L
a
m

b
d
a
 (

1
/P

a
)

Lambda vs A0 for OLS

1.962e-051.9618e-05 1.9622e-05

0e00

-5e-12

5e-12

Area (m2)

L
a
m

b
d
a
 (

1
/P

a
)

Lambda vs A0 for WLS U(Ap)

1.962e-051.9618e-05 1.9622e-05

0e00

-5e-12

5e-12

Area (m2)

L
a
m

b
d
a

 (
1

/P
a

)

Lambda vs A0 for WLS U2(Ap)

1.962e-051.9618e-05 1.9622e-05

0e00

-5e-12

5e-12

Area (m2)

L
a
m

b
d
a
 (

1
/P

a
)

Lambda vs A0 for GLS Ap

1.962e-051.9618e-05 1.9622e-05

0e00

-5e-12

5e-12

Area (m2)

L
a
m

b
d
a
 (

1
/P

a
)

Lambda vs A0 for WLS U(Ãp)

1.962e-051.9618e-05 1.9622e-05

0e00

-5e-12

5e-12

Area (m2)

L
a
m

b
d
a

 (
1

/P
a

)

Lambda vs A0 for WLS U2(Ãp)

1.962e-051.9618e-05 1.9622e-05

0e00

-5e-12

5e-12

Area (m2)

L
a
m

b
d
a
 (

1
/P

a
)

Lambda vs A0 for GLS Ãp

http://www.imeko.org/
https://www.scilab.org/