# Patent

## Title
Method for determining the flow rates of a first gas phase and at least one second liquid phase present in a polyphasic fluid

## Abstract

A method for determining a first flow rate of a gas phase and a second flow rate of at least one liquid phase present in a polyphasic fluid circulating in a conduit is presented. The method essentially comprises circulating the polyphasic fluid through the neck of a venturi delimited by the conduit. The method further comprises estimating the first flow rate and the second flow rate by using a measured pressure difference of the fluid through the neck and a property representative of the relative surface occupied by the gas phase relative to a total surface over a cross-section of the neck. The step of estimating the first flow rate and the second flow rate further comprises (a1) calculating an amount representative of the liquid phase present in the gas rich core relative to the total amount of liquid phase in the neck and calculating the first and second flow rate as a function of the amount representative of the liquid phase present in the core calculated during step (a1).

## Background
CROSS-REFERENCE TO RELATED APPLICATIONS
The present application is based on and claims priority to International Application No. PCT/FR2009/051806, filed 24 Sep. 2009, which is based on and claims priority to French Patent Application No. 0856466, filed 25 Sep. 2008. The entire contents of each are hereby incorporated herein by reference.
TECHNICAL FIELD
The present disclosure relates to a method for determining the first flow rate of a gas phase and the second flow rate of at least one liquid phase present in a polyphasic fluid circulating in a conduit.
BACKGROUND
During the exploitation of a well with view to producing hydrocarbons, it is known how to measure the flow rate of fluid extracted from the well in order to be able to control the amount and quality of the production.
However, the measurement of the fluid flow rate is complicated by the nature of the extracted fluid, which is generally polyphasic, with a gas phase and a liquid phase flowing at different rates.
It is therefore necessary for the operator of the well to determine the overall fluid flow rate flowing through the conduit and if possible, the individual volume flow rates of each phase flowing in the conduit.
SUMMARY OF THE DISCLOSURE
According to at least one aspect of the present disclosure, the method for determining the first flow rate of a gas phase and the second flow rate of at least one liquid phase present in a polyphasic fluid comprises the steps of:

    
    
        circulating the polyphasic fluid through a neck of a venturi delimited by a conduit, the polyphasic fluid substantially forming a gas rich core and a liquid-rich shell in the neck,
        measuring a pressure difference of the fluid through the neck;
        measuring a property representative of the relative surface occupied by the gas phase relative to a total surface over a cross-section of the neck;
        estimating the first flow rate and the second flow rate by using the measured pressure difference and the property representative of the relative surface occupied by the measured gas phase.
    
    


Such a method is intended to be applied in a polyphasic flow-meter. Such a flow-meter is notably used for characterizing the flow of a fluid extracted from a well made in the subsoil, like a well for producing hydrocarbons.
For this purpose, a polyphasic flow-meter of the type described in application US2006/0236779 can be used. Such a flow-meter comprises a venturi, a pressure tap through the venturi, and an apparatus for estimating the overall gas section present in a total section of the neck of the venturi. The flow-meter further comprises a computer for estimating the individual flow rates of the liquid phase and of the gas phase on the basis of measurements of pressure difference through the neck, of the measured or estimated relative surface, and on a calculation model. Such a flow-meter is efficient when the fluid contains a sufficiently large amount of liquid phase.
However, when the volume ratio of the gas phase over the total volume (GVF) is very high, notably when this ratio is greater than 90%, the accuracy of the flow-meter is degraded, in particular for the liquid flow. In this type of state of flow designated by the term of wet gas, the liquid phase is generally distributed in the form of an annular shell in contact with the wall of the conduit and the gas phase circulates at a greater rate in a central core delimited by the annular shell.
In one aspect of the present disclosure, a method is presented for determining the flow rates of a gas phase and of a liquid phase circulating through a polyphasic flow-meter which is very accurate, notably when the gas phase is a large majority.
For this purpose, the method may be characterized in that the step for estimating the first flow rate and the second flow rate comprises the following phases, or steps:
(a1) calculating an amount representative of the liquid phase present in the gas rich core relatively to the total amount of liquid phase in the neck depending on at least one amount representative of the flow of the polyphasic fluid in the neck and a first set of parameters depending on the geometry of the venturi and independent of the first flow rate and of the second flow rate,

(a2) calculating the first flow rate and the second flow rate as a function of the amount representative of the liquid phase present in the core calculated during step (a1).

The method according to the present disclosure may comprise one or more of the following features:

    
    
        the representative amount of the flow of the fluid circulating in the neck depends on the Weber number of the polyphasic fluid, and step (a1) comprises the calculation of the Weber number of the polyphasic fluid circulating in the neck;
        step (a1) may comprise the calculation of a corrected Weber number depending on the Weber number of the polyphasic fluid, and on the Reynolds number of the essentially liquid shell, the corrected Weber number forming the representative amount of the flow of the polyphasic fluid circulating in the neck;
        the method may comprise a preliminary step for calibrating the first set of parameters comprising the following phases:

(b1) successively putting into circulation through the neck, a plurality of polyphasic fluid samples having a first known gas phase flow rate and a second known liquid phase flow rate, each sample having at least a first flow rate and a second flow rate distinct from the other samples,

(b2) for each polyphasic fluid sample, calculating said or each amount representative of the flow of the polyphasic fluid in the neck, and the amount representative of the liquid phase present in the essentially gaseous core relatively to the total amount of liquid phase in the neck independently of said or each amount representative of the flow of the polyphasic fluid in the neck;

(b3) performing correlation between the representative amount of the liquid phase present in the essentially gaseous core relatively to the total amount of liquid phase in the neck and said or each amount representative of the flow of the polyphasic fluid in the neck in order to determine the parameters of the first set of parameters;

        the step for estimating the flow rates may comprise:

(c1) calculating a representative amount of the friction coefficient of the essentially liquid shell on the wall delimiting the venturi depending on a representative amount of the flow of the essentially liquid shell on the wall and of a second set of parameters depending on the geometry of the venturi and independent of the first flow rate and of the second flow rate;

        the representative amount of the flow of the essentially liquid shell is calculated as a function of the Reynolds number of the essentially liquid shell;
        step (b3) may comprise:

for each polyphasic fluid sample, the calculation of the representative amount of the flow of the essentially liquid shell on the wall and of the friction coefficient of the essentially liquid shell on the wall, independently of the representative amount of the flow of the essentially liquid shell on the wall,

the determination of the parameters of the second set;

the calculation of the representative amount of liquid phase in the essentially gaseous core as a function of the representative amount of the friction coefficient of the essentially liquid shell on the wall; and

the adjustment of the parameters of the first set of parameters for minimizing the error on the correlation between the representative amount of the liquid phase present in the essentially gaseous core relatively to the total amount of liquid phase in the neck and said or each representative amount of the flow of the polyphasic fluid in the neck;

        the step for estimating the flow rate comprises (c1) the calculation of a representative amount of the friction coefficient and the interface between the essentially gaseous core and the essentially liquid shell as a function of a representative amount of the flow of the essentially gaseous core in the neck, and of a second set of parameters depending on the geometry of the venturi and independent of the first flow rate and of the second flow rate;
        the calculation of the friction coefficient at the interface comprises the calculation of a simplified friction coefficient from the third set of parameters and from the Reynolds number of the essentially gaseous core in the neck and the calculation of a corrected friction coefficient representative of the irregularities of the interface between the liquid-rich shell and the essentially gaseous core on the basis of the simplified friction coefficient and of a fourth set of parameters;
        the step for estimating the first flow rate and the second flow rate comprises the calculation of a total mass flow rate through the venturi, depending on an overall density of polyphasic fluid, corrected as a function of the representative amount of liquid phase in the gas-rich core calculated in step (a1);
        the step for calculating the total mass flow rate comprises a step for calculating a constant proportional to the square root of the product of the corrected overall density and of an amount representative of the measured pressure difference, the constant being calculated as a function of the representative amount of the liquid phase present in the essentially gaseous core relatively to the total amount of liquid phase in the neck calculated during step (a1), and a representative amount of the friction coefficient between the wall and the essentially liquid shell; and
        the constant is further calculated depending on a representative amount of the compressibility of the essentially gaseous core circulating in the neck.
    
    





BRIEF DESCRIPTION OF THE DRAWINGS
To assist those of ordinary skill in the relevant art in making and using the subject matter hereof, reference is made to the appended drawings, which are not intended to be drawn to scale, and in which like reference numerals are intended to refer to similar elements for consistency. For purposes of clarity, not every component may be labelled in every drawing.
FIG. 1 is a sectional view along a median vertical plane of the venturi of the first polyphasic flow-meter in accordance with an exemplary embodiment disclosed herein.
FIG. 2 is a flow chart describing the main steps of the calibration phase in accordance with an exemplary embodiment disclosed herein.
FIG. 3 is a detailed flow chart of a particular step of the calibration phase described in the flow chart of FIG. 2.
FIG. 4 is a flow chart describing a phase for measuring individual volume flow rates of the gas phase and of the liquid phase during a measurement campaign in accordance with an exemplary embodiment disclosed herein.
FIG. 5 is a detailed flow chart of a particular step of the measurement phase described in the flow chart of FIG. 4.



DETAILED DESCRIPTION
In all the following, the terms “upstream” and “downstream” are meant relatively to the normal flow direction of a fluid in a conduit.
FIG. 1 illustrates a device 10 for measuring the flow rate of a fluid 12 flowing in a conduit 14 of a fluid exploitation installation such as a well for producing hydrocarbons.
The fluid 12 comprises a liquid phase and a gas phase. In this example, the volume ratio of the gas phase to the total volume of the fluid 12 flowing in the conduit 14, designated by the term of “Gas Volume Fraction” or GVF is advantageously greater than 90%. The fluid 12 thus comprises an essentially liquid annular shell 16 flowing at the contact of the wall delimiting the conduit 14 and an essentially gaseous core 18 flowing at the centre of the shell 16. These conditions are generally designated by the term of “annular flow”.
The conduit 14 for example extends vertically at the outlet of a well of the hydrocarbon exploitation installation (not shown). The fluid 12 flows in the conduit 14 along a vertical axis A-A′ facing the device 10. The essentially gaseous core 18 contains a portion of the liquid phase, notably as droplets 19 dispersed in the gas.
The conduit 14 delimits, in the vicinity of the measurement device 10, a venturi 20 comprising an upstream lower section 22 of inner diameter D, a downstream upper section 24 with an inner diameter substantially equal to D and, between the upstream section 22 and the downstream section 24, a venturi neck 26 with a diameter d smaller than D. The ratio β of the diameter d to the diameter D is for example comprised between 0.4 and 0.8.
The measurement device 10 comprises a sensor 28 for measuring the differential pressure Δp of the fluid between the upstream portion 22 and the neck 26, a sensor 30 for measuring the overall hold up of gas Γg and a computer 32 capable of estimating the total mass flow rate of the fluid ωt passing through the conduit 14 and the individual volume flow rates qg, ql of the gas phase and of the liquid phase, on the basis of a calculation model.
The sensor 30 for measuring the overall hold up of gas Γg, also designated by the term of “gas hold-up” comprises in this example a γ ray emission source 34 placed on one side of the conduit 14 in the neck 26 and a detector 36 for receiving the γ rays after their passing through the neck 26 into the fluid. The detector 36 is placed facing the source 34 on another side of the conduit 14. The γ rays emitted by the source transversely pass through the fluid 12 between the source 34 and the detector 36.
The sensor 30 provides the possibility of determining the overall linear gas hold up Γg corresponding to a homogeneous fluid flow by the equation (1) below:




  
    
      
        
          Γ
          g
        
        =
        
          
            1
            
              
                λ
                g
              
              -
              
                λ
                l
              
            
          
          ⁡
          
            [
            
              
                
                  1
                  d
                
                ⁢
                
                  ln
                  ⁡
                  
                    (
                    
                      
                        n
                        0
                      
                      n
                    
                    )
                  
                
              
              -
              
                λ
                l
              
            
            ]
          
        
      
    
    
      
        (
        1
        )
      
    
  




wherein n0 is the number of counts measured in an empty pipe in the absence of fluid, n is the number of counts received by the detector 36, λg is the linear attenuation of the pure gas, and λl is the linear attenuation of the pure liquid.

The computer 32 contains a model for calculating the first volume flow rate qg of the gas phase flowing through the conduit 14, the second volume flow rate ql of the liquid phase flowing through the conduit 14, on the basis of the measured pressure difference Δp, and of the overall linear gas hold up Γg, obtained by the sensors 28, 30.
As described in more detail below, the model is based on the calculation of the dispersed fraction ed of the liquid phase in the essentially gaseous core 18. This fraction ed is the ratio of the mass flow rate of the liquid which flows in the core 18 to the total mass flow rate of liquid which flows in the conduit 14.
According to an embodiment of the present disclosure, the ratio ed is calculated as a function of a representative amount of the flow of the polyphasic fluid 12, and to a first set of parameters p3, p4 which depend on the geometrical configuration of the venturi 20, independently of the liquid flow rates ql and of the gas flow rate qg passing through the conduit 14, by a first equation linking these quantities.
Advantageously, the ratio ed is calculated by the equation (2) below,

ed=p3×log(We′)+p4  (2)

wherein p3 and p4 are the parameters of the first set of parameters, and (We′) is a modified Weber number which depends on the Weber number (We) of the fluid 12 flowing in the neck 26, on the Reynolds number Ref of the essentially liquid shell 16 and on an adimensional density difference parameter εp, as defined by equation (3) and equation (4).

We′=We×εp−3/2×Ref1/2  (3)

ερ=(ρ1−ρg)/ρl  (4)

wherein ρl is the density of the liquid phase and ρg is the density of the gas phase.

The Weber number (We) is defined by the equation:




  
    
      
        We
        =
        
          
            2
            ⁢
            r
            ⁢
            
              
            
            ⁢
            
              ρ
              c
            
            ⁢
            
              j
              g
              2
            
          
          σ
        
      
    
    
      
        (
        
          4
          ⁢
          
            
          
          ⁢
          bis
        
        )
      
    
  




wherein r is the upstream radius of the venturi, σ is the surface tension between the gas and the liquid, ρc is the density of the core 18 and jg is the surface velocity of the gas.

The model further advantageously comprises the calculation of a wall friction coefficient cw which depends on a representative amount of the flow of the shell 16 in the neck 26, advantageously on the Reynolds number Ref of the shell 16, and on a second set of parameters p5, p6 which depend on the geometrical configuration of the venturi 20, independently of the liquid flow rate ql and gas flow rate qg passing through the conduit 14, by a second equation linking these quantities.
Advantageously, the wall friction coefficient Cw is calculated by the equation (5).

log(cw)=p5×log(Ref)+p6  (5)

The Reynolds number of the shell is given by the equation:




  
    
      
        
          Re
          f
        
        =
        
          
            2
            ⁢
            r
            ⁢
            
              
            
            ⁢
            
              ρ
              l
            
            ⁢
            
              u
              f
            
          
          
            η
            l
          
        
      
    
    
      
        (
        
          5
          ⁢
          
            
          
          ⁢
          bis
        
        )
      
    
  




wherein r is the radius of the neck of the venturi, uf is the average velocity of the film and ηl is the dynamic viscosity of the liquid.

The model further comprises the calculation of a gas friction coefficient cg which depends on an amount representative of the flow of the core 18, advantageously on the Reynolds number Rec of the gaseous core, and on a third set ph p2 of parameters which depend on the geometrical configuration of the venturi 20, independently of the liquid flow rate ql and gas flow rates qg passing through the conduit 14, by a third equation linking these quantities.
Advantageously, the coefficient Cg is calculated by the equation (6).

log(cg)=p1×log(Rec)+p2  (6)

The Reynolds number of the core is given by the equation:




  
    
      
        
          Re
          c
        
        =
        
          
            2
            ⁢
            
              
            
            ⁢
            
              ρ
              c
            
            ⁢
            
              
                h
                c
              
              ⁡
              
                (
                
                  
                    u
                    c
                  
                  -
                  
                    u
                    f
                  
                
                )
              
            
          
          
            η
            c
          
        
      
    
    
      
        (
        
          6
          ⁢
          
            
          
          ⁢
          bis
        
        )
      
    
  




wherein ρc is the density of the core, hc is the radius of the core upstream from the venturi, uc is the average velocity of the core, uf is the average velocity of the film, and ηc is the dynamic viscosity of the core.

In an advantageous embodiment, the model takes into account instabilities at the interface between the core 18 and the shell 16. In this embodiment, a friction coefficient Ci at the interface between the core 18 and the shell 16 is calculated by an equation linking it to the friction coefficient of the gaseous core cg through a fourth set of parameters w1, w2 which depend on the geometrical configuration of the venturi 20, independently of the liquid flow rates ql and gas flow rates qg passing through the conduit 14, through a fourth equation linking these quantities.
Advantageously this equation is defined below by:




  
    
      
        
          
            c
            1
          
          
            c
            2
          
        
        =
        
          1
          +
          
            
              w
              1
            
            ⁢
            
              
                
                  h
                  f
                
                d
              
            
            ⁢
            
              
                (
                
                  
                    u
                    g
                  
                  
                    u
                    gc
                  
                
                )
              
              
                w
                2
              
            
          
        
      
    
    
      
        (
        7
        )
      
    
  




wherein hf is the film height in the neck 26, d is the diameter of the neck 26, ug is the overall velocity of the gas in the core 18 and ugc is a critical gas velocity required for initiating wide irregular amplitude waves at the interface, as illustrated in the enlargement of FIG. 1.

The model present in the computer 32 is further based on writing a simplified sliding law which results from the equilibrium between the transfer of momentum at the interface between the core 18 and the shell 16, on the one hand, and the transfer of momentum at the interface between the shell 16 and the wall of the conduit 14, on the other hand, while neglecting the inertia and gravity terms.
This law may be written according to the equation:




  
    
      
        
          
            c
            w
          
          ⁢
          
            
          
          ⁢
          
            
              x
              f
              2
            
            
              
                α
                f
                2
              
              ⁢
              
                ρ
                l
              
            
          
        
        =
        
          
            c
            i
          
          ⁢
          
            
          
          ⁢
          
            
              ρ
              c
            
            
              
                α
                c
              
            
          
          ⁢
          
            
              (
              
                
                  
                    x
                    c
                  
                  
                    
                      ρ
                      c
                    
                    ⁢
                    
                      α
                      c
                    
                  
                
                -
                
                  
                    x
                    f
                  
                  
                    
                      ρ
                      l
                    
                    ⁢
                    
                      α
                      f
                    
                  
                
              
              )
            
            2
          
        
      
    
    
      
        (
        8
        )
      
    
  




wherein xf is the ratio of the liquid mass flow rate ωf in the shell 16 to the total mass flow rate ωt of fluid flowing in the conduit 14 and xc is the ratio of the mass flow rate ωc of liquid in the core 18 to the total mass flow rate ωt of fluid flowing in the conduit 14, αc is the effective hold up of the core and αf is the effective hold up of the film. The ratios of mass flow rates xf, xc themselves depend on the coefficient ed, as this will be seen below.

An exemplary method for determining the flow rate qg, ql according to the present disclosure will now be described with reference to FIGS. 2 to 5. This method initially comprises a calibration phase, illustrated by FIGS. 2 and 3 for determining the parameters p1 to p6, w1 and w2, for the given geometry of the conduit 14 on the basis of fluid samples with known flow rates. The method then comprises a phase for continuous measurement of the respective unknown individual flow rates ql and qg of the gas phase and of the liquid phase of a fluid 12 flowing in the conduit 14, illustrated by FIGS. 4 and 5.
In the calibration phase, a plurality of samples i of polyphasic fluids 12 preferably having a GVF ratio of more than 90% and a plurality of known liquid flow rates ql(i) and a plurality of known gas flow rates qg(i) are introduced into the conduit 14 so as to pass through the flow-meter device 10. For each known sample i, the respective density σl(i), σg(i) of the liquid phase and of the gas phase are determined experimentally as well as the dynamic viscosities ηl(i), ηg(i).
Next, as illustrated in step 50 in FIG. 2, a measurement of the pressure difference Δp(i) is carried out with the sensor 28 and a measurement of homogeneous gas hold up τg(i) is carried out with the sensor 30, for each sample corresponding to a known pair ql(i), qg(i).
Then, the parameters p1 to p6 and w1, w2 common to the whole of the samples, are determined by the computer 32. For this purpose, in step 52, the parameters of the third set p1, p2, and of the fourth set w1, w2 as defined in equations (5) and (6) and (7) are reset to a selected initial value. These values for example are p1=−1, p2=−0.5, w1=0, w2=1. The values of ci and cw are also reset to a given value for example equal to 0.005.
Next, in step 54, an iterative loop for optimizing the parameters p5, p6 of the second set is performed, as illustrated by FIG. 3. At each iteration, the loop 54 comprises in step 56, a step for calculating the ratio ed(i) by a sliding law for each pair of flow rates (ql(i); qg(i))
As stated earlier, this sliding law is based on the equilibrium between the mass transfers between the interface and the wall according to equation (8).
In equation (8), the ratios xf and xc are defined by the equations:

xf=ωf/ωt  (9)

and

xc=ωc/ωt  (10)

αc is the effective hold up of the gas in the core 18 as defined by the equation:





  
    
      
        
          α
          c
        
        =
        
          
            
              α
              g
            
            ⁡
            
              (
              
                1
                +
                
                  
                    e
                    d
                  
                  ⁢
                  
                    
                  
                  ⁢
                  
                    
                      q
                      l
                    
                    
                      q
                      g
                    
                  
                
              
              )
            
          
          =
          
            
              α
              g
            
            ⁡
            
              (
              
                1
                +
                
                  
                    e
                    d
                  
                  
                    Q
                    g
                  
                
              
              )
            
          
        
      
    
    
      
        (
        11
        )
      
    
  




and ρc is given by the equation:





  
    
      
        
          ρ
          c
        
        =
        
          
            
              ρ
              g
            
            (
            
              1
              +
              
                
                  e
                  d
                
                ⁢
                
                  
                
                ⁢
                
                  
                    1
                    -
                    
                      
                        e
                        d
                      
                      ⁢
                      
                        
                          ρ
                          g
                        
                        
                          ρ
                          l
                        
                      
                    
                  
                  
                    
                      
                        ω
                        g
                      
                      
                        ω
                        l
                      
                    
                    +
                    
                      
                        e
                        d
                      
                      ⁢
                      
                        
                      
                      ⁢
                      
                        
                          ρ
                          g
                        
                        
                          ρ
                          l
                        
                      
                    
                  
                
              
            
            )
          
          =
          
            
              ρ
              g
            
            ⁡
            
              (
              
                1
                +
                
                  
                    e
                    d
                  
                  ⁢
                  
                    
                  
                  ⁢
                  
                    
                      1
                      -
                      
                        
                          e
                          d
                        
                        ⁢
                        
                          R
                          g
                        
                      
                    
                    
                      
                        X
                        g
                      
                      +
                      
                        
                          e
                          d
                        
                        ⁢
                        
                          R
                          g
                        
                      
                    
                  
                
              
              )
            
          
        
      
    
    
      
        (
        
          11
          ⁢
          
            
          
          ⁢
          bis
        
        )
      
    
  




wherein αg is calculated from the gas hold up Γg in a homogeneous fluid by the equation (12).





  
    
      
        
          α
          g
        
        =
        
          
            
              Γ
              g
              2
            
            ·
            
              (
              
                1
                +
                
                  
                    e
                    d
                  
                  ⁡
                  
                    (
                    
                      
                        q
                        l
                      
                      
                        
                          q
                          g
                        
                        ⁢
                        
                          
                        
                      
                    
                    )
                  
                
              
              )
            
          
          =
          
            
              Γ
              g
              2
            
            ⁡
            
              (
              
                1
                +
                
                  
                    e
                    d
                  
                  
                    Q
                    g
                  
                
              
              )
            
          
        
      
    
    
      
        (
        12
        )
      
    
  




On this basis, equation (8) may be written according to the equation




  
    
      
        
          
            c
            w
          
          
            c
            i
          
        
        =
        
          
            
              
                1
                +
                
                  1
                  /
                  
                    A
                    g
                  
                
              
              
                1
                +
                
                  
                    e
                    d
                  
                  /
                  
                    Q
                    g
                  
                
              
            
          
          ·
          
            
              
                e
                d
              
              +
              
                
                  R
                  g
                
                ⁢
                
                  Q
                  g
                
              
            
            
              
                e
                d
              
              +
              
                Q
                g
              
            
          
          ·
          
            
              (
              
                
                  
                    
                      
                        Q
                        g
                      
                      /
                      
                        A
                        g
                      
                    
                    -
                    
                      e
                      d
                    
                  
                  
                    1
                    -
                    
                      e
                      d
                    
                  
                
                -
                1
              
              )
            
            2
          
        
      
    
    
      
        (
        
          12
          ⁢
          
            
          
          ⁢
          bis
        
        )
      
    
  




wherein Qg is the ratio of the flow rates defined by the equation

Qg=qg/ql  (12ter)

Ag is the ratio of the hold ups as defined by the equations

Ag=αg/αl  (13)





  
    
      
        
          A
          g
        
        =
        
          
            
              Γ
              g
              2
            
            ⁡
            
              (
              
                1
                +
                
                  
                    e
                    d
                  
                  /
                  
                    Q
                    g
                  
                
              
              )
            
          
          
            1
            ⁢
            
              
                –Γ
                g
                2
              
              ⁡
              
                (
                
                  1
                  +
                  
                    
                      e
                      d
                    
                    /
                    
                      Q
                      g
                    
                  
                
                )
              
            
          
        
      
    
    
      
        (
        
          13
          ⁢
          
            
          
          ⁢
          bis
        
        )
      
    
  




Rg is the ratio defined by the equation:

Rg=ρg/ρl  (14)

and Ng is the ratio defined by equation (15),

Ng=ηg/ηl  (15)

wherein ηg is the dynamic viscosity of the gas and ηl is the dynamic viscosity of the liquid.

Once ed(i) is calculated for each pair of known values qg(i), ql(i) the amount xf and the amount αf are calculated in step 58.
For this purpose, the amount xf is determined by the equation:




  
    
      
        
          x
          f
        
        =
        
          
            
              x
              l
            
            ×
            
              e
              d
            
          
          =
          
            
              
                
                  ω
                  l
                
                
                  ω
                  t
                
              
              ·
              
                e
                d
              
            
            =
            
              
                1
                -
                
                  e
                  d
                
              
              
                
                  1
                  +
                  
                    
                      R
                      g
                    
                    ·
                    
                      Q
                      g
                    
                  
                
                ⁢
                
                  
                
              
            
          
        
      
    
    
      
        (
        16
        )
      
    
  




and the amount αf is determined by the equation:

αf=1−αc  (16bis)

wherein αc is calculated by equation (11) and αg is calculated by equation (12).

This having been done, for each sample I, in step 60, the quantity cw(i) is calculated from the equation:




  
    
      
        
          ω
          
            t
            ⁢
            
              
            
          
        
        =
        
          C
          ·
          ɛ
          ·
          
            a
            1
          
          ·
          
            
              
                
                  2
                  ⁢
                  
                    ρ
                    ~
                  
                
                
                  1
                  -
                  
                    β
                    4
                  
                
              
              ⁢
              Δ
              ⁢
              
                
              
              ⁢
              
                p
                ′
              
            
          
        
      
    
    
      
        (
        17
        )
      
    
  




for each sample i, which links the total flow rate ωt of fluid in the conduit 14, to a corrected pressure difference Δp′, to an estimated global density {tilde over (σ)}, via a proportionality constant C·ε at the section of the neck a1.

This equation results from the integration of the Navier Stokes equations over the length of the venturi.
Equation (18) gives the estimated global density {tilde over (ρ)} as a function of xf and of xc,




  
    
      
        
          1
          
            ρ
            ~
          
        
        =
        
          
            
              x
              f
              2
            
            
              
                ρ
                l
              
              ⁢
              
                α
                f
              
            
          
          +
          
            
              x
              c
              2
            
            
              
                
                  ρ
                  c
                
                ⁢
                
                  α
                  c
                
              
              ⁢
              
                
              
            
          
        
      
    
    
      
        (
        18
        )
      
    
  




Equation (19) gives the dynamic pressure difference Δp′ as a function of the measured pressure difference Δp and of a corrective parameter in order to take gravity into account.

Δp′=Δp+(αfρl+αcρc)gΔz  (19)

Equation (20) gives the proportionality constant in which C is the discharge coefficient and ε is the compressibility factor.




  
    
      
        
          
            (
            
              1
              
                C
                ·
                ɛ
              
            
            )
          
          2
        
        =
        
          1
          +
          
            
              1
              
                1
                -
                
                  β
                  4
                
              
            
            ⁢
            
              (
              
                
                  
                    c
                    ω
                  
                  ⁢
                  
                    κ
                    ω
                  
                  ⁢
                  
                    
                  
                  ⁢
                  
                    
                      x
                      f
                      2
                    
                    
                      α
                      f
                      2
                    
                  
                  ⁢
                  
                    
                      ρ
                      ~
                    
                    
                      ρ
                      l
                    
                  
                
                +
                
                  
                    ξ
                    g
                  
                  ⁢
                  
                    
                      ρ
                      ~
                    
                    
                      ρ
                      
                        g
                        ,
                        0
                      
                    
                  
                
              
              )
            
          
        
      
    
    
      
        (
        20
        )
      
    
  




In these equations, Δz is the height between two points of the pressure tap 28 and κw is the surface of the wall divided by the volume of the venturi and ρg,0 is the density of the gas upstream from the venturi.
The term ξg is defined by equation (21) in order to take into account the compressibility of the gas during its passing through the neck.




  
    
      
        
          
            ξ
            g
          
          
            
              a
              1
              2
            
            ⁢
            
              ρ
              
                g
                ,
                0
              
            
          
        
        =
        
          
            β
            4
          
          ⁢
          
            
              δ
              P
            
            κ
          
          ⁢
          
            
              x
              c
              2
            
            
              α
              c
            
          
          ⁢
          
            
              X
              g
            
            
              
                X
                g
              
              +
              
                e
                d
              
            
          
          ⁢
          
            
              ∫
              1
              
                1
                /
                
                  β
                  4
                
              
            
            ⁢
            
              
                
                  -
                  t
                
                ·
                
                  d
                  t
                
              
              
                
                  [
                  
                    1
                    -
                    
                      
                        δ
                        P
                      
                      ⁡
                      
                        (
                        
                          1
                          -
                          t
                        
                        )
                      
                    
                  
                  ]
                
                
                  1
                  +
                  
                    1
                    /
                    κ
                  
                
              
            
          
        
      
    
    
      
        (
        21
        )
      
    
  




wherein:

Xg=ωg/ωl  (22)

δp=(1−p1/p2)/(1/β4−1)  (23)

The term κ is the polytropic coefficient calculated with the gas and liquid mass flow rates ωg, ωl, the specific capacities of the gas and of the liquid cv,g, cv,l, and the isentropic exponent γ from equation (24), and the terms p1 and p0 are the respective pressures at the neck and upstream from the neck.




  
    
      
        κ
        =
        
          
            
              
                ω
                g
              
              ⁢
              
                c
                
                  v
                  ,
                  g
                
              
              ⁢
              γ
            
            +
            
              
                ω
                l
              
              ⁢
              
                c
                
                  v
                  ,
                  l
                
              
            
          
          
            
              
                ω
                g
              
              ⁢
              
                c
                
                  v
                  ,
                  g
                
              
            
            +
            
              
                ω
                l
              
              ⁢
              
                c
                
                  v
                  ,
                  l
                
              
            
          
        
      
    
    
      
        (
        24
        )
      
    
  




Once cw (i) is calculated for each pair ql(i), qc(i), the Reynolds numbers Rec(i) of the core 18 and Ref(i) of the shell 16 are calculated in step 62.
The Reynolds number Rec(i) of the core is calculated as a function of ed by the equation:




  
    
      
        
          Re
          c
        
        =
        
          
            2
            
              π
              ⁢
              
                
              
              ⁢
              
                r
                0
              
            
          
          ⁢
          
            
              ρ
              g
            
            
              η
              g
            
          
          ⁢
          
            q
            g
          
          ⁢
          
            
          
          ⁢
          
            
              
                
                  
                    Q
                    g
                  
                  ⁢
                  
                    R
                    g
                  
                
                +
                
                  e
                  d
                
              
              
                
                  
                    Q
                    g
                  
                  ⁢
                  
                    R
                    g
                  
                
                +
                
                  
                    e
                    d
                  
                  /
                  
                    N
                    g
                  
                
              
            
            ·
            
              
                
                  
                    Q
                    g
                  
                  ⁢
                  
                    R
                    g
                  
                
                +
                
                  e
                  d
                
              
              
                
                  
                    Q
                    
                      g
                      ⁢
                      
                        
                      
                    
                  
                  ⁢
                  
                    R
                    g
                  
                
                +
                
                  
                    e
                    d
                  
                  ⁢
                  
                    R
                    g
                  
                
              
            
            ·
            
              1
              
                
                  Γ
                  g
                  2
                
                ⁡
                
                  (
                  
                    1
                    +
                    
                      
                        e
                        d
                      
                      /
                      
                        Q
                        g
                      
                    
                  
                  )
                
              
            
          
        
      
    
    
      
        (
        25
        )
      
    
  




The Reynolds number Ref(i) of the shell 16 is calculated for each sample pair i as a function of ed by the equation:




  
    
      
        
          Re
          f
        
        =
        
          
            2
            
              π
              ⁢
              
                
              
              ⁢
              
                r
                0
              
            
          
          ⁢
          
            
              ρ
              l
            
            
              
                η
                l
              
              ⁢
              
                
              
            
          
          ⁢
          
            
              q
              l
            
            ⁡
            
              (
              
                1
                -
                
                  e
                  d
                
              
              )
            
          
          ⁢
          
            
              
                1
                /
                
                  A
                  g
                
              
              -
              
                
                  e
                  d
                
                /
                
                  Q
                  g
                
              
            
            
              1
              +
              
                1
                /
                
                  A
                  g
                
              
            
          
        
      
    
    
      
        (
        26
        )
      
    
  




wherein r0 is the upstream radius of the neck.

A plurality of pairs (cw(i); Ref(i)) associating the friction coefficient of the shell on the wall of the conduit 14 with the Reynolds number of the shell are obtained.
In step 64, a linear regression is made on the equation:

log(cw)=p5×log(Ref)+p6  (27)

for calculating the coefficient for p5 and p6 on the basis of the pairs (cw(i); Ref(i)) associated with each pair of flow rates.

Next, in step 66, the coefficient ci(i) is calculated by the equation (12 bis) as a function of Rec(i), ed(i) Qg(i), et Γg(i). For this purpose Rec(i) is calculated as a function of ed(i) and of Qg(i) by equation (25). Next, the coefficient cg(i) is calculated on the basis of the correlation of equation (6) as a function of Rec(i) and of coefficients p1 and p2.
Next, with equation (11), the coefficient αc is calculated. The coefficient




  
    
      h
      f
    
    d
  




is then determined by the equation:





  
    
      
        
          
            
              h
              f
            
            d
          
        
        =
        
          
            
              1
              -
              
                
                  α
                  c
                
              
            
            2
          
        
      
    
    
      
        (
        
          27
          ⁢
          
            
          
          ⁢
          bis
        
        )
      
    
  




The coefficient ci(i) is then calculated on the basis of the coefficients w1 and w2 and of equation (7).
In step 68, the differences |Δp5| and |Δp6| between the coefficients p5 and p6 obtained in step 64 during this iteration of the loop 54 and those obtained during the preceding iteration of the loop are calculated.
If at least one of the differences |Δp5| and |Δp6| is greater than a given convergence coefficient, for example 10−6, a new iteration is carried out by returning to step 56. If the value of each of these differences is less than the given convergence coefficient, the loop 54 is completed and step 70 is applied. In step 70, the error εw committed on the coefficient cw during the linear regression made in step 64 is estimated.
This error εw is for example calculated by equation (28).




  
    
      
        
          ɛ
          w
        
        =
        
          
            1
            N
          
          ⁢
          
            
              
                
                  ∑
                  
                    i
                    =
                    1
                  
                  N
                
                ⁢
                
                  log
                  ⁡
                  
                    (
                    
                      C
                      
                        w
                        ,
                        i
                      
                    
                    )
                  
                
              
              -
              
                (
                
                  
                    
                      p
                      5
                    
                    ×
                    
                      log
                      ⁡
                      
                        (
                        
                          Re
                          
                            f
                            ,
                            i
                          
                        
                        )
                      
                    
                  
                  +
                  
                    p
                    6
                  
                
                )
              
            
          
        
      
    
    
      
        (
        28
        )
      
    
  




In step 72, an optimization test of this error εw is carried out. If the error εw is always greater than a given optimized value, the coefficients p1, p2, w1 and w2 are modified in step 74, for example by a descent along the gradient.
A new iteration of the loop 54 is then carried out for calculating new coefficients p5, p6, by using the coefficients p5, p6 obtained during the preceding iteration for initializing the loop in step 56. When the error εw is less then the given optimized value, the optimization loop 54 of the coefficients p1, p2, w1, w2 is then stopped. The coefficients ed(i) are then recovered for each flow rate pair ql(i), qg(i) and the modified Weber number, We′(i) as calculated by equation (3) is calculated on the base of the Weber number determined by means of the equation:




  
    
      
        We
        =
        
          
            2
            
              
                π
                2
              
              ⁢
              
                r
                0
                3
              
            
          
          ⁢
          
            
              ρ
              g
            
            σ
          
          ⁢
          
            
              q
              g
              2
            
            ·
            
              
                
                  
                    Q
                    g
                  
                  ⁢
                  
                    R
                    g
                  
                
                +
                
                  e
                  d
                
              
              
                
                  
                    Q
                    g
                  
                  ⁢
                  
                    R
                    g
                  
                
                +
                
                  
                    e
                    d
                  
                  ⁢
                  
                    R
                    g
                  
                
              
            
          
        
      
    
    
      
        (
        29
        )
      
    
  




The coefficients p3 and p4 are then calculated by linear regression in step 76.
A complete set of parameters p1 to p6 and w1, w2 is therefore obtained during the calibration phase. This allows the quantities ed, cg, cw, ci to be calculated as a function of these parameters and of characteristic quantities of the fluid flow through the conduit during a measurement phase, according to the equations (2), (5), (6) and (7), as this will be detailed below.
Once the calibration phase is completed, a measurement phase is carried out with a fluid 12 of unknown flow rate circulating in the conduit 14. This measurement may be carried out periodically at intervals for example comprised between 5 minutes and 15 minutes.
As earlier described, the pressure difference Δp is measured by the sensor 28 and the overall gas hold up Γg is calculated by the equation (1) with the measurement made by the sensor 30 in step 80.
Next, the parameters ed, cw and ci are initialized by having them assume a value given in step 82. This value is for example equal to 0.5 for ed, 0.05 for cw, and 0.08 for ci respectively.
Next, a loop 84 of iterations is carried out for determining the respective flow rates ql and qg. This loop 84 begins by a step for calculating the ratio Qg, as defined by equation (12ter) calculating it on the basis of the equation (12bis), in step 86.
Next, in step 88, an iteration loop is carried out for calculating the total mass flow rate ωt. This loop 88 is described in FIG. 5. It begins with a phase for calculating ωt with equation (17), taken as a combination with equations (18) to (20) (step 90) by using the values of ed, cw, and ci.
In step 92, the liquid volume flow rate ql is calculated by the equation:

ql=ωt/(1+1/(RgQg))  (30)

and in step 94, the Reynolds number of the film Ref is calculated by equation (26).

Next, in step 96, the coefficient c, is calculated by equation (5) by using the parameters p5 and p6 determined during the calibration phase.
Next, in step 98, a convergence test is carried out on the difference |Δωt| between the value of ωt calculated in step 90 during this iteration of the loop 84 and the value ωt calculated during the preceding iteration of the loop 84.
If this difference |Δωt| is greater than a convergence value given for example as equal to 10−6, a new iteration is accomplished on the loop 84 by returning to step 90 and by using the value of cw calculated in step 96. However, if this difference |Δωt| is less than the convergence value, the loop 84 is completed and the values of ωi and cw are then extracted in step 100.
Next, with reference to FIG. 4 in step 102, the corrected Weber number We′ is calculated by equation (29), taken in combination with equations (3), (4) and (26).
Next, the Reynolds number Rec of the core 18 is calculated by equation (25), taken in combination with the equation

qg=ωt/(ρ1(1+RgQg))  (31)

In step 104, the coefficients ed and ci are respectively calculated by the correlations defined by equation (2) and by the combination of equations (6) and (7) in which the parameters p1 to p6 and w1, w2 are those calculated during the calibration phase. In step 106, the equations (30), (31) are used for calculating the flow rates qg and ql. In step 108, a test is carried out on the respective differences |Δql| and |Δqg| between the values qg, ql calculated in step 106 during this iteration of the loop 84 and the respective values qg, ql calculated during the preceding iteration of the loop 84.
If at least one of these differences |Δql| and |Δqg| is greater than a determined convergence coefficient, for example equal to 10−6, a new iteration of the loop 84 is carried out by returning to step 86 and by using the new obtained values of ed, cw and ci. If these differences |Δql| and |Δqg| are smaller than the given convergence coefficient, the loop 84 is stopped and the coefficient ql and qg are extracted so as to be for example displayed by the computer 32, with the value of the overall mass flow rate ωt.
In a first alternative, the interface between the essentially liquid shell 16 and the essentially gaseous core 18 is considered as smooth. In this case the friction coefficient at the interface ci is equal to the gas friction coefficient cg and the parameters w1 and w2 are equal to 0 during all the steps of the method.
In another alternative, the gas compressibility ξg in equation (20) is considered as 0 during the whole method, so that equations (21) to (24) are not used. The errors obtained on the total mass flow rate ωt and on the individual gas and liquid flow rates qg, ql, for a fluid circulating at a pressure of more than 25 bars in the conduit 14 through the measurement device 10 are summarized in Table 1 below.










TABLE 1









Total mass flow rate (wt)
±1.6%



Gas volume flow rate (qg)
±2.1%












Liquid volume flow rate (ql)
90%, 
90%, 
90%, 




<GVF 
<GVF 
<GVF 




<97%
<99%
<100%




±3.2%
±7%
±10 m3/day










As illustrated by this table, using a model estimating the ratio ed of the liquid mass circulating in the core 18 as a function of a number We′ representative of the flow of the fluid 12 through the conduit 14, it is possible to obtain excellent accuracy on the estimated value of the gas volume flow rate qg and on the liquid volume flow rate ql, even for high GVF values and close to 100%.
This estimation of the quantity ed, made in combination with the estimation of the friction coefficient at the interface ci between the core 18 and the shell 16 and of the friction coefficient cw between the shell 16 and the wall of the conduit 14 also contributes to improving the accuracy of the measurement.
By accurately defining the total mass flow rate ωt by equations (18) to (20), it is also possible to better take into account physical phenomena occurring in the conduit 14 for improving the accuracy of the measurement.

## Claims

1. A method for determining a first flow rate (qg) of a gas phase and a second flow rate (q) of at least one liquid phase present in a polyphasic fluid circulating in a conduit, the method comprising:
circulating the polyphasic fluid through a neck of a venturi delimited by a conduit, the polyphasic fluid forming a gas-rich core and a liquid-rich shell in the neck,
measuring a pressure difference (Δp) of the fluid through the neck;
measuring a property (Γg) representative of the relative surface occupied by the gas phase relative to a total surface over a cross-section of the neck;
estimating the first flow rate (qg) and the second flow rate (ql) by using the measured pressure difference (Δp) and the property (Γg) representative of the relative surface occupied by the measured gas phase;
wherein estimating the first flow rate (qg) and the second flow rate (ql) further comprises:
(a1) calculating an amount (ed) of the liquid phase present in the gas-rich core relative to the total amount of liquid phase in the neck depending on at least one amount (We′) representative of the flow of the polyphasic fluid in the neck and a first set of parameters (p3, p4) depending on the geometry of the venturi and independent of the first flow rate (qg) and of the second flow rate (q),
(a2) calculating the first flow rate (qg) and the second flow rate (ql) as a function of the amount (ed) representative of the liquid phase present in the core calculated during (a1).

2. The method according to claim 1, wherein the amount (We′) representative of the flow of the fluid circulating in the neck is a function of the Weber number (We) of the polyphasic fluid, and (a1) comprises the calculation of the Weber number (We) of the polyphasic fluid circulating in the neck.

3. The method according to claim 2, wherein (a1) further comprises the calculation of a corrected Weber number (We′) as a function of the Weber number (We) of the polyphasic fluid, and of the Reynolds number (Ref) of the liquid shell, the corrected Weber number (We′) forming the amount representative of the flow of the polyphasic fluid circulating in the neck.

4. The method according to claim 1, wherein the method further comprises a preliminary step for calibrating the first set of parameters (p3, p4), wherein the preliminary step comprises:
(b1) successively circulating through the neck of the venturi a plurality of polyphasic fluid samples having a first known gas phase flow rate (qg(i)) and a second known liquid phase flow rate (ql(i)), each sample having at least one first flow rate or a second flow rate distinct from the other samples;
(b2) for each polyphasic fluid sample, calculating each amount (We′) representative of the flow of the polyphasic fluid in the neck, and of the amount (ed) representative of the liquid phase present in the gaseous core relative to the total amount of liquid phase in the neck independently of each amount (We′) representative of the flow of the polyphasic fluid in the neck; and
(b3) performing correlation between the amount (ed) representative of the liquid phase present in the gaseous core relatively to the total amount of liquid phase in the neck and each amount (We′) representative of the flow of the polyphasic fluid in the neck for determining the first set of parameters (p3, p4).

5. The method according to claim 1, wherein estimating the first flow rate (qg) and the second flow rate (ql) further comprises:
(c1) calculating an amount representative of the friction coefficient (cw) of the liquid shell on the wall delimiting the venturi as a function of an amount (Ref) representative of the flow of the liquid shell on the wall and of a second set of parameters (p5, p6) depending on the geometry of the venturi and independent of the first flow rate and of the second flow rate.

6. The method according to claim 5, wherein the amount (Ref) representative of the flow of the liquid shell is calculated as a function of the Reynolds number of the liquid shell.

7. The method according to claim 4, wherein for each polyphasic fluid sample, (b3) comprises:
calculating an amount (Ref) representative of the flow of the liquid shell on the wall and of the friction coefficient (cw) of the liquid shell on the wall, independently of the amount (Ref) representative of the flow of the liquid shell on the wall,
determining a second set of parameters (p5, p6) depending on the geometry of the venturi and independent of the first flow rate and of the second flow rate;
calculating the amount (ed) representative of the liquid phase in the gaseous core as a function of the amount representative of the friction coefficient (cw) of the liquid shell on the wall;
adjusting the first set parameters (p3, p4) for minimizing the error in the correlation between the amount (ed) representative of the liquid phase present in the gaseous core relative to the total amount of liquid phase in the neck and each amount (We′) representative of the flow of the polyphasic fluid in the neck.

8. The method according to claim 1, wherein estimating the first flow rate (qg) and the second flow rate (ql) further comprises:
(c1) calculating an amount representative of the friction coefficient (ci) at the interface between the gaseous core and the liquid shell as a function of an amount (Rec) representative of the flow of the gaseous core in the neck, and of a third set of parameters (p1, p2) depending on the geometry of the venturi and independent of the first flow rate (qg) and of the second flow rate (ql).

9. The method according to claim 8, wherein the calculation of the friction coefficient at the interface (ci) comprises the calculation of a simplified friction coefficient (cg) from the third set of parameters (pi, p2) and of the Reynolds number (Rec) of the gaseous core in the neck and the calculation of a corrected friction coefficient (ci) representative of the irregularities of the interface between the liquid-rich shell and the gaseous core on the basis of the simplified friction coefficient (cg) and a fourth set of parameters (w1, w2).

10. The method according to claim 1, wherein the step for estimating the first flow rate (qg) and the second flow rate (ql) comprises the calculation of a total mass flow rate (ωt) through the venturi as a function of an overall polyphasic fluid density ({tilde over (σ)}) corrected according to the amount (ed) representative of liquid phase in the gas-rich core calculated in (a1).

11. The method according to claim 10, wherein calculating the total mass flow rate (ωt) comprises calculating a constant C·ε proportional to the square root of the product of the corrected overall density ({tilde over (σ)}) and of an amount representative of the measured pressure difference (Δp), the constant C·ε being calculated as a function of the amount (ed) representative of the liquid phase present in the gaseous core relative to the total amount of liquid phase in the neck calculated during (a1), and an amount (cw) representative of the friction coefficient between the wall and the liquid shell.

12. The method according to claim 11, wherein the constant C·ε is further calculated as a function of an amount representative of the compressibility (ξg) of the gaseous core circulating in the neck.