International Journal of Multiphase Flow

Slug flow is a common flow pattern in pipelines that is often accompanied by undesirable effects like vibrations, pressure loss, and corrosion. Since these effects correlate with slug frequency, various attempts to predict this parameter by empirical or semi-empirical methods have been undertaken in the past. However, significant mismatches between these predictions can be observed. In this work, different slug frequency calculation methods have been applied to simulation data to investigate the sensitivity of threshold parameters that are often used in slug detection algorithms. The findings reveal that the detection of slugs from liquid holdup data is highly sensitive to these thresholds. Aeration of the liquid phase causes the gas–liquid interface to be less distinct and requires an adaption of the thresholds to the degree of aeration. In contrast, slug detection algorithms based on frequency analysis are robust to small deviations of the liquid level but fail to properly discriminate between slugs and waves. Our investigations show that slug frequency strongly depends on the method chosen for the determination of the liquid level. We propose new approaches that are less susceptible to aeration and approximate the liquid level very close to the authors’ human judgment.


Introduction
The slug flow pattern is one of the most commonly observed twophase flow patterns in industrial applications (Pineda-Pérez et al., 2018).For horizontal flow, it is characterized by alternating blocks of aerated liquid (so-called slugs) and gas bubbles flowing above liquid films (Al-Safran, 2009).These liquid slugs can lead to an increased pressure drop, vibrations, as well as increased corrosion and erosion of the pipe (Hill and Wood, 1994;Al-Safran, 2009;Marcano et al., 1998).Hence, they can cause severe problems in industrial operations and lead to large uncertainties in multiphase flow metering (Olbrich et al., 2021a).Because these unfavorable effects strongly correlate to slug frequency (Al-Safran, 2009), this parameter is of special interest.
Due to the inherent unsteady and random behavior of slug flow, it is considered as one of the most complex flow patterns (Marcano et al., 1998), and slug frequency is difficult to predict.Various attempts with empirical and semi-empirical approaches can be found in the literature, see Gregory and Scott (1969), Heywood and Richardson (1979), Hill and Wood (1990), Jepson and Taylor (1993), Al-Safran (2009), Schulkes (2011).However, different slug frequency correlations predict significantly different slug frequencies for the same flow conditions.Even the predicted trends can differ: For the test cases considered in this paper, some correlations predict a slug frequency that decreases As established above, these criteria for identifying slugs are crucial for the comparison of slug frequencies determined by different investigations.However, detailed information about the applied criteria are scarce in the literature: Gregory and Scott (1969), Jepson and Taylor (1993), as well as Woods and Hanratty (1999), Woods et al. (2006) identified slugs from pressure pulses but did not mention an exact pressure difference or threshold that had to be surpassed by the slugs to be counted.The latter two (Woods and Hanratty, 1999;Woods et al., 2006) acknowledge the described difficulties and argue that pressure pulses are the preferable method for slug detection because the slug bodies get aerated and it becomes difficult to distinguish waves and slugs by liquid level or liquid holdup at high superficial gas velocities.Marcano et al. (1998), Gokcal et al. (2009), Al-Safran (2009), and Zabaras (2000) utilized capacitance/conductance sensors to determine slug frequencies from liquid holdup, but also without clarifying the exact criteria for a signal to be interpreted as a slug.Moreover, two of these (Al-Safran, 2009;Zabaras, 2000) gathered slug frequency data from the literature, and it is not clear, whether they considered different measurement and calculation methods among their respective data sets.Another instrumentation was utilized by Heywood and Richardson (1979), who measured the liquid holdup by means of the -ray absorption method and then performed a power spectral density (PSD) analysis to determine the most dominant frequency, which was considered to be the slug frequency.El-Oun (1990) also utilized rays to determine slug frequency, but the calculation method was not described.In contrast, Hernandez-Perez et al. (2010) declared that a liquid holdup larger than 0.7 is a common criteria for a slug structure and cited Nydal (1991) as a reference.Instead of a fixed liquid holdup threshold, Zhao et al. (2013) as well as Baba et al. (2017) applied a variable threshold, namely the middle between the highest and lowest measured liquid holdup in a time series.More recently, Soedarmo et al. (2019) and Soto-Cortes et al. (2021) developed an automatic threshold selection procedure based on liquid holdup probability density functions for pseudo-slug and slug flow, respectively.They applied a two-threshold approach to identify slug front and tail, which then enabled the calculation of slug length, and hence, slug liquid holdup.
Despite the known difficulties, only little attention has been paid to the exact characteristics of common measurement and calculation methods and their impact on the slug frequency.For high viscosity fluids, Zhao et al. (2013) compared the slug frequencies derived by a thresholding and a PSD method and found good agreement.However, the comparison was limited to low superficial gas velocities, where aeration is low and the shape of the plugs/slugs is more distinct compared to higher superficial gas velocities.This leads to the questions how different slug detection algorithms as well as the parameters used in these methods affect the resulting slug frequency and how thresholding methods compare to frequency analyses by means of PSD at higher superficial gas velocities.
To shed more light on these questions, we simulated horizontal slug flow at seven different superficial gas velocities with the opensource software package OpenFOAM v1812.The numerical simulations provide data for the whole 3D domain.Hence, they give insight into areas that are hardly accessible by experiments.The spatially and temporally resolved liquid volume fraction data enables us to apply and compare different evaluation methods and criteria for the considered test cases.Slug frequency is calculated in two steps: First, the liquid level is estimated from the liquid volume fraction data.Four different approaches are investigated herein.Second, the slug frequency is calculated from the liquid level time series.Here, two fundamentally different methods, namely frequency analysis via PSD and liquid level thresholds, are considered.For the latter one, we conduct a thorough sensitivity analysis of the thresholds that need to be defined for these approaches.
The present work shows how the different calculation methods influence the slug frequency.This clarifies for which cases it is necessary to know the exact criteria of how a slug frequency was determined,  Mandhane et al. (1974).
and for which part of the plug/slug flow regime the slug frequency is robust and insensitive to the utilized methods.The investigations help to judge slug frequency predictions from the literature that have been determined with different methods and to choose appropriate measurement and calculation methods for slug frequencies in future research.

Simulation methodology
This section describes the general simulation set-up, the utilized schemes, the boundary and initial conditions, as well as a mesh convergence study.We consider the flow of Paraflex oil and nitrogen through a horizontal pipe of inner diameter  = 0.097 m and length  = 305 .The fluid properties prescribed in the numerical simulations are summarized in Table 1.Seven different test cases are considered.All test cases have the same superficial liquid velocity,   = 1.87 m∕s, whereas the superficial gas velocity is different for each case:   [m∕s] ∈ {0.35, 0.62, 1.12, 1.87, 2.5, 3.0, 4.0}.According to the flow pattern map of Mandhane et al. (1974), they all lie within the plug or slug flow regime, where the transition from plug to slug flow occurs at around   = 0.75 m∕s for the considered liquid superficial velocity   = 1.87 m∕s, see Fig. 1.
For the numerical simulation, the open-source software package OpenFOAM v1812 was used.From this package, the two-phase solver interFoam was chosen, which is suitable for the transient simulation of two incompressible, isothermal immiscible fluids.It is based on the volume of fluid (VOF) method (Hirt and Nichols, 1981), which is a numerical technique for tracking and locating the interface between two fluids.Since the flow is turbulent, an unsteady Reynolds-averaged Navier-Stokes (URANS) approach was applied.Menter's shear stress transport (SST) turbulence model (Menter, 1993) was employed to close the system of equations.
For the spatial discretization of the pipe geometry, a block-structured mesh was used, see Fig. 2. The refinement level needed for the numerical simulations was determined by means of a grid study.In this study, simulations with a superficial gas velocity of 0.63 m/s were performed for eight different meshes.The results were compared to each other with respect to slug frequency since this is the parameter of interest in this paper.The mesh was refined separately in radial  and longitudinal direction.In radial direction, two refinement levels were considered:   = 46 ∕ and   = 54 ∕, where   denotes the number of cells per diameter.In longitudinal direction, four different refinement levels were compared:   [1 ∕m] ∈ {85, 100, 120, 150} with   denoting the number of cells per meter in longitudinal direction.To save computational costs, only a pipe length of  = 165  instead of the 305 D was considered in the grid study, which is still sufficiently long for the formation of slugs.Furthermore, only half of the pipe geometry was simulated.For this, the computational domain was virtually cut along the vertical symmetry plane in longitudinal direction (-plane).
A symmetryPlane boundary condition was applied to the new face, which was treated like a slip boundary condition (no friction).
Fig. 3 shows the comparison of slug frequency over distance to the inlet for the different mesh resolutions.One can see that the mesh refinement in radial direction has hardly any influence on the resulting slug frequency (compare the solid and dotted lines for each color).The mesh refinement in longitudinal direction, on the other hand, leads to increasing slug frequencies, at least for the first two refinement levels.Nevertheless, grid convergence can be observed for   = 120 ∕m and   = 150 ∕m.Altogether, it can be concluded that a mesh resolution with   = 46 cells in radial and   = 120 cells per meter in longitudinal direction is sufficient.
For four selected meshes (1.  = 46 ∕,   = 85 ∕m; 2.   = 46 ∕,   = 100 ∕m; 3.   = 54 ∕,   = 85 ∕m; and 4.   = 54 ∕,   = 120 ∕m), additional simulations for the full geometry (i.e., the whole cylinder) were performed to investigate the effect of simulating only half of the pipe geometry and using symmetry boundary conditions.The investigations showed that, for the coarser meshes, differences in the resulting slug frequencies can be observed.Using the full mesh leads to higher slug frequencies.Hence, it has a similar effect as a grid refinement in longitudinal direction.For the most finest mesh considered in this part of the study (  = 54 ∕,   = 120 ∕m), these differences become much smaller.In this case, the mean relative error of the slug frequency determined by the half mesh with respect to the full mesh is less than 2.9 %.A brief investigation of required computation time revealed that the shape of the mesh (half or full cylinder) and the resolution in radial direction had a significant impact on computational resources.In contrast, increasing the resolution in longitudinal direction was less demanding.Considering the results of the mesh convergence study and the evaluation of computational resources, a mesh of half-cylindrical shape with a resolution of   = 46∕ and   = 14∕ ≈ 144 ∕m was used for the numerical simulations.The final calculations were run on a high performance cluster using 24 processors for each simulation.Depending on the prescribed superficial gas velocities, a simulation of 70 s took between one and four weeks (higher superficial gas velocities were computationally more demanding).
For the time discretization, the implicit Euler scheme was used.The time step size was adjusted automatically by limiting the Courant number to 0.5.The total simulation time was 70 s.The schemes used for the spatial discretization are summarized in Table 2. Regarding initial conditions, at the start of the simulation, the pipe contains no liquid and an initial velocity of 0 m/s is defined.On the pipe walls, no-slip boundary conditions are applied.The inlet cross section was bisected horizontally as shown in Fig. 2. The gas enters the pipe through the upper part, the liquid through the lower part of the inlet.To enhance slug initiation, the inlet perturbation method by Schmelter et al. (2021a) is used.In this approach, the secondary velocity components (parallel to the inlet face) are changed over time.They are chosen randomly in direction and amplitude, but only up to 20 % of the main component's amplitude.Since only the secondary velocity components are disturbed, the total volume flux is not influenced.This approach enhances the formation of slugs in the pipe, which is beneficial for computational costs because the length of the domain can be reduced.Furthermore, in Schmelter et al. (2021a) it is shown that the mean slug frequency at downstream locations converges for different perturbation amplitudes.Hence, it can be concluded that the inlet perturbations do not affect the slug frequency in developed slug flow.

Extraction of liquid level time series from simulation results
To calculate the slug frequency   , the liquid level is analyzed.According to Andritsos and Hanratty (1987), it is defined as the vertical position of the gas-liquid interface relative to the inner pipe diameter in the vertical centerline of the pipe cross section.This means that ''liquid level'' refers to the dimensionless liquid level ℎ  = ℎ  ∕, where ℎ  is the dimensioned liquid level, i.e., the height of the liquid-gas interface in the middle of the pipe.However, determining the liquid level is not trivial.In this section, we describe several methods to approximate the liquid level.

Approximation of liquid level based on area-average of liquid volume fraction
The liquid holdup   at a fixed longitudinal position  in the pipe is calculated as the area-average of the liquid volume fraction (, , , ) for all time points  as (1) Here,  denotes the cross-sectional area of the pipe: where  is the inner diameter of the pipe.
Under the assumption of a plain, horizontal interface between liquid and gas, the liquid holdup   is equivalent to the area of a circular segment and can be used to approximate the liquid level ℎ area  by solving the following equation iteratively: where  = ∕2 is the radius of the cross section.The liquid holdup is a parameter that is often used in measurements to determine slug frequencies (e.g., in Marcano et al. (1998), Heywood and Richardson (1979), Gokcal et al. (2009), Zhao et al. (2013), Baba et al. (2017)).By means of Eq. ( 2) such measurements of   can easily be converted into a corresponding liquid level ℎ area  .The conversion into ℎ area  allows a better comparison with the liquid levels introduced in the following.

Approximation of liquid level based on line-average of liquid volume fraction
Instead of averaging the liquid volume fraction over the whole cross-sectional area , in the following approach, it is only averaged over the vertical centerline of the pipe cross section: To calculate ℎ line  numerically,  = 39 probes are located equidistantly with a distance of ∕( + 1) = ∕40 to each other and to the walls, see Fig. 4.
The liquid level is then approximated by: where  1 and   are the outmost probe locations.These probes represent the distance between the probe and the wall of ∕( + 1) plus half of the distance to the next probe of ∕(2( + 1)).Thus, they cover 1.5 times more distance than the inner probes at locations   ,  = 2, … ,  − 1.
In a separated flow with a level interface, this would yield the liquid level with a tolerance of the distance between the probes of ∕( + 1).However, due to aeration of the liquid phase or liquid droplets in the gas phase, the liquid level determined by this method can be lower or higher than the actual liquid level.Bubbles and droplets that intersect with one of the  probes have a larger impact on the liquid level determined by this method compared to the previous method based on the liquid holdup.
To account for bubbles included in the liquid phase, which are typical for slug flow, the following correction algorithm is proposed.In a first iteration, ℎ line  is calculated.In a second step, an auxiliary liquid volume fraction  impr is calculated by setting all liquid volume fraction values of probes located below ℎ line  to one.Then, the average of the auxiliary liquid volume fraction is calculated.Altogether, ℎ impr  is determined similar to (4) using  impr instead of : where  impr (  , , ),  = 1, … , , is given by In a truly separated flow, no bubbles occur and all probes below ℎ line , are shown.The time intervals were chosen to demonstrate the representative characteristics of each method.In an ideally separated flow, the interface and, therefore, the liquid level could be determined unambiguously at the point of transition between alpha values of zero and one.This is nearly the case for   = 0.35 m/s, where the phases are clearly separated.Hence, all three methods align well.For   = 1.87 m/s, on the other hand, significant differences between the liquid level approximations can be observed at 14.6 s, where a gas bubble is enclosed in the liquid phase.(Note that the shape of the bubble is distorted and squashed by this visualization).In this case, ℎ line

Approximation of liquid level based on binary methods
The idea of this approach is to convert the continuous liquid volume fraction field  into a field of binary values (either 0 or 1) by means of a threshold  bin .If  >  bin , the fluid is counted as liquid.Otherwise, it is considered as gas.The liquid level ℎ bin  at time  is then defined as the highest -position occupied by liquid: Using the previously defined probes at locations   ,  = 1, … , , it can be approximated as follows:  0.35 m/s, 1.87 m/s, and 4.0 m/s, respectively, where the time interval is again chosen to demonstrate the characteristics of the compared methods.For low aeration with a clear interface (  = 0.35 m∕s), the different thresholds  bin = 0.5 and  bin = 0.8 yield nearly the same results.However, when the superficial gas velocity is increased, individual droplets cause the ℎ bin  to jump for  bin = 0.5.Because the droplets usually have a low liquid volume fraction value, the graph for  bin = 0.8 remains smooth and reproduces the shape of the liquid structure better.When the superficial gas velocity is increased further, as shown in the plot with   = 4 m/s, the number of spikes caused by individual droplets increases for  bin = 0.5.In contrast, the graph for  bin = 0.8 exhibits downward spikes because of aeration.However, while an individual droplet in the gas phase is sufficient to cause an upward spike, an individual bubble in the liquid phase cannot cause a downward spike (because of the maximum used in the definition of ℎ bin  ).Therefore, the downward spikes that are typical for high thresholds  bin occur more rarely and are less pronounced.

Calculation of slug frequency from liquid level time series
The slug frequency is the number of slugs that passes a specific point along the pipeline over a certain period of time (Al-Safran, 2009).In this work, we considered a time interval of 50 s (the initial 20 s of the simulation were omitted) and evaluated the slug frequency at position  = 300  downstream from the inlet.For the detection of slugs from the liquid level time series, two different approaches are used.The first one is based on defining thresholds for slug detection, the second one uses frequency analysis.

Thresholds for slug detection
Soto-Cortes et al. ( 2021) and Soedarmo et al. (2019) applied a two threshold method to identify the slug front and slug tail, respectively.This allows to calculate the slug length and suppresses double counting of the same slug structure.For the present investigation, a similar methodology with two thresholds is applied.However, a simpler approach is used.This could be done because the thresholds discussed herein are solely evaluated to determine slug frequency but not slug length.When the liquid level exceeds the upper threshold  up a slug is counted.The lower threshold  low makes sure that small fluctuations of the liquid level do not lead to double counting of slugs.
Choice of upper threshold.The upper threshold  up defines the required minimum liquid level of a structure to be identified as a slug.A typical characteristic of slugs is that the liquid phase occupies the whole cross section.Thus, for ideal slugs, the liquid level should increase to one.Hence, for the methods that compensate aeration, like ℎ bin an upper threshold  up =0.95 (i.e., a value close to 1) is chosen.On the other hand, a much lower value for  up is needed for methods that do not take aeration into account.To determine an appropriate value for  up for these methods, the aeration of the slugs is estimated by means of slug holdup prediction methods from the literature.A review of several methods is presented in Appendix A. If not stated otherwise, we use the slug holdup prediction by Gomez et al. (2000).Eq. ( 2) is utilized to convert the liquid holdup of the slugs into a corresponding liquid level.As mentioned early on, not all researchers describe the exact method how they calculated the slug frequencies.However, using liquid level or liquid holdup thresholds as explained here appears to be a popular method, see Hernandez-Perez et al. ( 2010), Baba et al. (2017), Zhao et al. (2013), Fabre et al. (1995).
Choice of lower threshold.The lower threshold  low can be understood as a threshold to identify the tail of a slug body.Bubbles or an uneven liquid surface may cause the liquid level of a slug body to fluctuate and cross  up several times.The lower threshold assures that these fluctuations are still accounted to the same slug body until the liquid level falls below  low .It is chosen as the value of the average liquid level.As shown below in Section 5.1, this value is a good choice for slug detection.For the 2.5 m/s case, the liquid level does not exceed  up between 21.0 s and 21.5 s.However, the liquid volume fraction plot reveals that the liquid phase is aerated and reaches the top of the pipe.This demonstrates how aeration affects the calculated liquid level and why  up has to be adapted to the superficial gas velocity or slug holdup.However, each slug is unique and may have a different shape and level of aeration.Therefore, a certain  up might be too high to detect all slugs, but too low to also count waves as slugs at the same time.
Comparing the liquid volume fraction plots on the right, one observes how aeration increases with higher superficial gas velocities.In the plots of the liquid level time series (left picture), on the other hand,  a change from relatively smooth curves to curves with a lot of peaks and oscillations can be observed, indicating the transition from plug to slug flow at a superficial gas velocity of around 1.12 m/s.This value is in a similar range as the value predicted by the flow pattern map of Mandhane et al. (1974), see Fig. 1.

Frequency analysis
Another approach to determine slug frequency from liquid level time series is frequency analysis.It is, for example, applied in Hernandez-Perez et al. ( 2010), Knotek et al. (2016), Schmelter et al. (2021a).One advantage of frequency analysis is that no thresholds need to be chosen in advance.Furthermore, frequency analysis methods (like a Fourier transform) are easy to apply.On the other hand, the most dominant frequencies determined by these methods are not necessarily the frequencies of the slugs, but only reflect the dynamics of the interface between the different phases (Schmelter et al., 2021a;Knotek et al., 2021).
In this paper, Welch's PSD estimate (function pwelch in Matlab) was applied to the liquid level time series at position  = 300 .
The function pwelch calculates the one-sided PSD estimate using Welch's segment averaging estimator, see The MathWorks, Inc. ( 2019) for details.Due to additional smoothing, it is more robust than a pure Fourier transform.Table 3 summarizes the parameters used for the computation of pwelch.In the simulation, data were saved every 0.0025 s leading to a sampling frequency of 400 Hz.Hence, the segment length of the filter corresponds to a time interval of 1.5 s.The number of points used in the discrete Fourier transform (DFT) corresponds to a time interval of 15 s.

Results and discussion
In the following, we discuss the sensitivity of slug frequency to different slug detection algorithms as well as to the parameters used in these methods.In Section 5.1, the influence of thresholds on the slug frequency is investigated.Then, in Section 5.2, the effect of the prediction of the liquid holdup in slugs is discussed.Finally, a comparison of the different methods is presented in Section 5.3.

Sensitivity of slug frequency to thresholds
As established above, a lot of slug detection algorithms use thresholds to decide whether a structure is counted as slug or not.In the following, the sensitivity of the resulting slug frequency with respect to these thresholds is presented.
In Fig. 8, the slug frequency determined from ℎ area  is plotted against  up .The solid lines represent the slug frequencies determined with  low equal to the average liquid level, whereas the dotted lines show the slug frequency determined with  low equal to  up , which is equivalent to not utilizing a lower threshold.The circles on each solid curve mark the thresholds based on the slug holdup estimation by Gomez et al. (2000).This estimation is described in detail in Appendix A. If a lower threshold  bin is used, one observes hardly any sensitivity to  up for the cases in the plug flow regime, i.e., superficial gas velocities of 0.35 m/s and 0.62 m/s.The corresponding lines are close to constant for  up < 0.95.On the other hand, for cases with higher superficial gas velocity, the resulting slug frequency depends strongly on the choice of  up .Thus, the choice of  up is crucial for the determination of the slug frequency.
Because  up has to be larger than  low and  low is different for each case, the solid curves for each case are defined on different intervals.The comparison of the dotted lines with the corresponding solid lines exhibits the effect of utilizing a lower threshold.The lower threshold always reduces the determined slug frequency, which is particularly significant for low superficial gas velocities.This is because the lower threshold prevents fluctuations of the liquid level of a slug body to be counted as slugs, see top picture in Fig. 7. Hence, a lower threshold is necessary to distinguish between several slug structures on one hand and fluctuations of the liquid level within one slug structure on the other hand.For ℎ line  and ℎ impr  (not shown here), similar observations can be made as for ℎ area  .For ℎ impr , the key difference is that all slug frequencies are slightly increased, and, therefore, all curves are shifted to the right.
In summary, the slug frequency is highly sensitive to  up .Furthermore,  low is efficient to avoid the double counting of slugs due to fluctuations in the liquid level within one slug structure.
Fig. 9 shows the slug frequency determined from ℎ area  plotted against  low to visualize the sensitivity of the slug frequency to  low .The plotted intervals are not identical for each curve because  low has to be lower than  up , and all curves are evaluated for different values of  up corresponding to the slug holdup predicted by Gomez et al. (2000).The average liquid level of each case, which has been chosen as ''default value'' for the lower threshold in this paper, is marked with a dot on each curve.The chosen thresholds generally lie within an interval of almost constant slug frequency, especially for cases with either high or low superficial gas velocity.The largest influence of  low on the resulting slug frequency is observed for cases with medium superficial gas velocities   [m/s] ∈ {1.12, 1.87, 2.5}.However, the influence is still small compared to  up .This shows that small deviations of the chosen  low have a low impact on the determined slug frequency.Again, similar trends are observed for ℎ line , where ℎ impr  generally yields higher slug frequencies.
Table 4 shows how the slug frequency changes when the thresholds  low and  up are increased or decreased by 0.05. low has the largest impact on the case with a superficial gas velocity of 1.12 m/s, where the slug frequency is reduced by 17.9 %, when  low is reduced.In average decreasing or increasing  low by 0.05, changes the slug frequency by −5.2 % and +3.9 %, respectively.
In contrast, a deviation of the chosen  up by 0.05 increases or decreases the slug frequency by at least 33.9 % for cases with superficial gas velocities of 1.12 m/s or higher.In average, decreasing  up by 0.05 increases the slug frequency by 92.1 %, whereas reducing  up reduces the slug frequency by 46.3 %.This shows that the choice of  up has a large impact on the resulting slug frequency, at least for the cases in the slug flow regime (i.e., superficial gas velocity of 1.12 m/s and higher), whereas the slug frequency is much less sensitive to  low .Fig. 10 shows the sensitivity of the slug frequency with respect to the binary threshold  bin used in ℎ bin  .As a reminder, the slug frequencies from ℎ bin  are all calculated with an upper threshold  up = 0.95 since aeration is taken into account by ℎ bin  .For each case, the slug frequency exhibits a decreasing trend for increasing  bin .Setting  bin = 0.5 seems intuitive to separate the gas from the liquid phase.However, as previously shown in Fig. 6, droplets in the gas phase may cause significant fluctuations of ℎ bin  for low values of  bin .For the cases with relatively low superficial gas velocities (  [m/s] ∈ {0.35, 0.62, 1.12}), slug frequency is almost equal for  bin = 0.5 and  bin = 0.8.Most cases in the slug flow regime show a (light) dependency on  bin in this range, presumably because droplets are less miscounted as slugs for higher thresholds  bin .Therefore, a threshold  bin = 0.8 is chosen, and the corresponding slug frequencies are marked with dots in the figure .In Fig. 10, one can also observe that the slug frequency decreases with increasing superficial gas velocity and that the   -over- bin curves do not cross each other.Only the 1.12 m/s case, which is close to the     plug-to-slug flow transition as observed in Fig. 7, is an exception to both observations.For superficial gas velocities of 2.5 m/s or larger, Fig. 8 revealed that the slug frequency drops to 0 /s for  up = 0.95.This indicates that ℎ area  never exceeds values of 0.95 for these cases.In contrast, Fig. 10 shows slug frequencies greater than 0 for  bin = 0.95 for all superficial gas velocities.This shows that the liquid phase reaches almost the top of the pipe, but aeration reduces the liquid level predicted by ℎ area  .

Sensitivity of slug frequency to liquid holdup prediction
As aforementioned, an estimation of the slug holdup is required to calculate the slug frequency from ℎ area  and ℎ line  , and we exemplarily utilized the prediction by Gomez et al. (2000) in the previously presented investigations.However, a variety of different slug holdup estimations can be found in the literature.Some of them are briefly presented in Appendix A, namely Gregory et al. (1978), Malnes (1982), Ferschneider (1983), Andreussi et al. (1993), Marcano et al. (1998), Gomez et al. (2000), andAbdul-Majeed (2000).Note that, the selection of methods is based on the two review papers by Pereyra et al. (2012) and Ibarra et al. (2019).To investigate how the choice of the slug holdup correlation affects the resulting slug frequency, we calculated the slug frequency from ℎ area  multiple times, each time with a different  up according to each slug holdup correlation.Fig. 11 shows the seven different predictions of the slug frequency for the seven different slug holdup correlations spanning a large range of slug frequencies.Andreussi et al. (1993) predict a slug liquid holdup of 0.998 for   = 0.35 m∕s (see Appendix A).A value of this magnitude can only be achieved by ideal plugs/slugs so that it results in a slug frequency of 0 /s in the simulation data.All other slug liquid holdup predictions result in a similar slug frequency of about 1.5 (slugs) per second at   = 0.35 m∕s.For   = 0.62 m∕s, the slug frequency decreases to values between about 1 and 1.35 (slugs) per second for all slug holdup prediction methods.If the superficial gas velocity is further increased, however, the slug frequencies do not show the same trend anymore, but diverge and, hence, span a remarkable range.This is in line with the observations made in Fig. 8 that the slug frequency is less sensitive in the plug flow region, but aeration and lower slug liquid holdups make slug frequency more sensitive to the upper threshold for higher superficial gas velocities in the slug flow region.This is particularly troublesome, when the upper threshold is based on other sensitive quantities like the slug liquid holdup.
Summed up, Fig. 11 reaffirms that the thresholding methods are highly sensitive to the chosen upper threshold (that is based on slug liquid holdup estimations), especially in the slug flow region.This means that two slug frequency predictions based on the same data set can differ completely if thresholds based on different holdup estimations are chosen.

Comparison of different slug frequency calculation methods
In previous sections, four different approximations of the liquid level (Sections 3.1-3.3)as well as two methods to calculate the slug frequency from the liquid level were described (Sections 4.1-4.2).Therefore, eight combinations of liquid level approximations and slug detection algorithms are possible.In Fig. 12, these eight slug frequency predictions are plotted over the superficial gas velocity.The four liquid level approximations (ℎ bin  , ℎ impr  , ℎ area  , ℎ line  ) are discerned by different colors, while the slug detection algorithms (thresholding and PSD) are plotted in solid and dashed lines, respectively.For the thresholding  algorithms utilizing ℎ area  and ℎ line  , a variable upper threshold based on the liquid holdup by Gomez et al. (2000) was applied, as described earlier in Section 4.1.Conversely, a fixed upper threshold of 0.95 was applied for ℎ bin  and ℎ area

𝑙
because those two liquid level approximations compensate for aeration.
Furthermore, several slug frequency predictions from the literature (Al-Safran, 2009;Gregory and Scott, 1969;Heywood and Richardson, 1979;Jepson and Taylor, 1993;Schulkes, 2011;Zabaras, 2000) are plotted by gray dotted lines with markers.For more details about these prediction methods, the reader is referred to Appendix B, where the markers are assigned to their designated references.
Fig. 12 allows several observations: first, all eight determined slug frequencies exhibit a decreasing trend over superficial gas velocity and align well with each other in the plug flow regime, where aeration is low and high liquid levels are obtained over a sustained period of time.In contrast, the slug frequencies determined by PSD and thresholding behave differently in the slug flow regime, and the differences increase with superficial gas velocity.At superficial gas velocities between   = 2.5 m∕s and   = 4 m∕s, two groups can be recognized: The slug frequencies determined by the PSD show a generally increasing trend of slug frequency over superficial gas velocity, whereas the thresholding method results in a falling trend.
This first observation partly matches with the results reported by Zhao et al. (2013).An important difference is that Zhao et al. (2013) reported a good match between PSD and thresholding up to superficial gas velocities of 2 m/s, we observed large discrepancies between the methods already superficial gas velocities of 1.12 m/s.However, Zhao et al. (2013) performed their comparison for superficial liquid velocities of 0.1 m/s, which, according to the flow pattern map shown in Fig. 1, shifts the plug-to-slug flow transition to higher superficial gas velocities.Hence, their results are in good agreement with the observations made in this study.
Second, the four PSD curves match each other closely, which implies that the considered methods of liquid level approximation is less relevant for the slug frequency and that the simplest approximation may suffice in this case.In contrast, the thresholding method obtains a wide spread of slug frequencies for the different liquid level approximations.
And third, among the thresholding methods, the slug frequencies based on ℎ bin  and ℎ impr  match closely, which are the two approximations that compensate for aeration and aim to detect the highest point of the liquid phase.In contrast, ℎ area  and ℎ line  , which are based on averages of the liquid volume fraction, both yield different results compared to the other methods.These deviations can be explained by the inconsistent predictions of the liquid holdup from the literature, see Fig. 11.
The thresholding methods aim to identify actual slug structures, while the PSD is only able to identify dominant structures, but it is not clear whether these dominant structures are indeed slugs or only waves.It is a well-known drawback of methods based on frequency analysis that they cannot clearly distinguish between waves and slugs (Schmelter et al., 2020;Knotek et al., 2021;Schmelter et al., 2021b).The slug frequency based on thresholding of ℎ bin  is not expected to underestimate the slug frequency, because a single probe located at  ≥ 0.95  measuring a liquid volume fraction of ≥ 0.8 (see Section 3.3) is sufficient to trigger a slug being counted.Therefore, we assume that the thresholding methods for ℎ area  and ℎ line  , and especially the PSD method (for all four liquid level approximations) overestimate the slug frequency in the slug flow regime.This shows that the height of the liquid structure cannot reliably be determined either by PSD or liquid holdup.Hence, the methods that aim to detect the top of the liquid phase, i.e., ℎ bin  and ℎ impr  , are less susceptible to aeration and perform best.
Another observation in Fig. 12 is made by comparing the calculated slug frequencies with predictions from literature.Among the predicted slug frequencies, the values and trends of slug frequencies differ remarkably.Monotonic increasing, monotonic decreasing, as well as non-monotonic slug frequency predictions can be found in the literature.For example, the prediction by Heywood and Richardson (1979) exhibits a non-monotonic trend and is based on slug frequencies derived by PSD, which is in line with our results.The slug frequency predictions also reveal a high uncertainty throughout the shown interval, where slug frequencies range from 0.25 to 1.25 slugs per second, for instance at the slug-to-plug flow transition.For medium superficial gas velocities (0.65 m/s ≤   ≤ 2.5 m/s), our results lie within the range of predicted slug frequencies.Of note, many slug frequency predictions are based on mixed data published or provided by several individual research teams.One can readily see that the average of our eight calculated slug frequencies would fit in the predictions very well.

Conclusions
Difficulties in slug frequency prediction have been observed in the literature.Possible reasons are the random and unsteady nature of slug flow, the large variety of parameters that influence slug frequency as well as inaccurate liquid holdup measurements, e.g., by tomography (Olbrich et al., 2021b).However, also the approaches used to determine the holdup or liquid level in the pipe as well as the methods used to calculate the slug frequency have an influence on the resulting slug frequency.
To assure that the applied method is feasible to determine the slug frequency, researchers often verify their approach by comparison with F. Webner et al. visual observations.However, human perception is also prone to error, especially for fast movements.The slugs passing by create droplets or a liquid film on the inner side of the glass viewing section, which further hinders observation.Van Hout et al. (2003) stated that even image processing with the help of video cameras is of limited use because the gas-liquid interface is hard to detect for highly aerated slugs.Soedarmo et al. (2019) andSoto-Cortes et al. (2021) proposed objective methods to avoid arbitrary parameter selection in slug detection.
However, besides Zhao et al. (2013), who compared slug frequencies determined with thresholding against those determined by PSD, little attention has been paid to the influence and bias of different calculation methods.Such a comparison is provided in this paper.
To apply different slug frequency calculation methods to the same data set, seven different superficial gas velocities at constant superficial liquid velocity (  = 1.87 m∕s) and constant diameter ( = 0.097 m) were simulated with OpenFOAM v1812.Four approximations of the liquid level were applied, including ℎ area  , which is equivalent to evaluating the liquid holdup as often done in experiments.Then, two different approaches were followed to calculate the slug frequency from the liquid level approximations, namely the thresholding and frequency analysis approach.Therefore, a total of eight combinations of liquid level approximation and slug frequency calculation methods were investigated and following observations were made: 1.In the plug flow regime, the plugs show a clear defined shape, the plug liquid holdup is close to one by definition (Yadigaroglu and Hewitt, 2017), and all methods yield similar slug frequencies.Deviations of the slug frequency among the methods are much smaller than in the slug flow regime.2. The PSD is not sensitive to the method of liquid level determination.However, only the dominant structures are analyzed and, hence, it is not feasible to distinguish between slugs and waves.3.In the slug flow region, slugs get aerated with increasing superficial gas velocity.Because each method is affected by aeration differently, the methods result in different slug frequencies.Moreover, when approximating the liquid level by ℎ line  or ℎ area  , an estimation of the aeration or slug liquid holdup is required, which significantly increases the uncertainty of the prediction.4. The presented binary methods overcome the problem of adapting the threshold to the degree of aeration.However, they have the disadvantage that droplets in the gas phase might be miscounted as slugs.5.The proposed ℎ impr  also accounts for aeration and represents the liquid level better than the other methods.However, it requires detailed knowledge about the liquid volume fraction field data.
Without knowing the exact level of aeration, the height of the waves/slugs cannot be reliably determined from the liquid holdup.Therefore, the methods aiming to directly detect the top of the liquid phase, i.e., ℎ bin  and ℎ impr  , performed better in slug detection, especially in distinguishing slugs from waves.Since PSD and liquid holdup both fail to reliably distinguish aerated slugs from waves, methods that can detect whether the liquid structure breaches to the top of the pipe, like ℎ impr  , are recommended to identify aerated slug structures.Altogether, the paper reveals that the choice of the evaluation method has a significant impact on the resulting slug frequency, especially in the slug flow region.This insight is necessary to judge and compare slug frequencies determined by different research teams utilizing different calculation methods.Furthermore, it helps to choose calculation methods for future research and emphasizes the necessity of standardized slug counting criteria for highly aerated slugs.
To gain further insight into the influence of slug counting criteria on slug frequency, more data with different superficial liquid velocities, diameters, and fluid properties is required.While numerical simulations provide a detailed view of the whole flow field for individual cases, experimental setups are recommended to generate results for a large variety of flow conditions.Furthermore, experiments could validate the behavior of the presented slug detection algorithms for measured data.

Declaration of competing interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.(1978), who investigated 2.54 cm and 5.12 cm internal diameter pipes.Jepson and Taylor (1993) 2000) compared existing slug frequency prediction methods to a data set of 399 data points.194 of those data points were collected from published literature that covered internal diameters ranging from 2.54 cm to 20.32 cm.Slug frequency data of various fluid properties were included in the data set and the inclination angle was varied between 0 • and 11 • .Zabaras (2000) adapted the approach from Gregory and Scott (1969) to the collected data and found the following correlation: where   and   are the actual liquid velocity and gas velocity, respectively.To calculate the actual velocities of the fluids, the occupied cross section of each phase is required.The procedure by Taitel and Dukler (1976) is utilized to calculate the liquid level and the gas and liquid velocity.Schulkes (2011) investigated the influence of several parameters on slug frequency.Published data with a total of about 1200 data points were used.He found that fluid properties like density and viscosity can be neglected for turbulent flow, but should be considered for laminar flow.Furthermore, he found no significant correlation between pressure and slug frequency.Schulkes (2011)  Even though one would expect an equal or similar slug frequency at same flow conditions, a large range of slug frequencies is covered by the different prediction methods.Several different trends can be observed: The methods by Jepson and Taylor (1993) and Al-Safran (2009) exhibit a monotonically increasing trend, whereas the curve by Schulkes (2011) decreases monotonically.The prediction methods by Heywood and Richardson (1979) and Zabaras (2000) are inspired by Gregory and Scott (1969).Hence, these three methods show a similar trend, with a local minimum between 2 m/s and 3 m/s.Further analyses, differences and similarities of different slug frequency prediction methods from the literature can be found in Webner (2021).

F
.Webner et al.

Fig. 4 .
Fig. 4. Illustration of the probes used for the calculation of the liquid level ℎ line  .

𝑙
automatically attain the value of one.In this case, ℎ impr  = ℎ line  .In case of slug flow with a lot of aeration, on the other hand, ℎ impr  > ℎ line  .Fig. 5 shows three plots of the liquid volume fraction over time at the  probes for three different superficial gas velocities:   [m∕s] ∈ {0.35, 1.87, 4.0}.Each dot represents the liquid volume fraction () of the respective time and -position.The higher the liquid volume fraction, the darker the dot.A black dot represents the liquid phase ( = 1) and a white one the gas phase ( = 0).Furthermore, the approximations of the liquid level determined by the three different methods, ℎ area  , ℎ line  , and ℎ impr and ℎ area  both underestimate the liquid level, whereas ℎ impr  is able to predict it correctly.At around 14.52 s, droplets in the air lead to a slight overestimation of the liquid level by ℎ line  .This overestimation is further increased by ℎ impr  .For   = 4.0 m/s, the aeration increases and the differences between the averaging methods compared to ℎ impr  become more apparent: The averaging methods never capture the actual height of the wave structure, whereas ℎ impr  reproduces the shape of the wave significantly better.At around 30.69 s, ℎ impr  exhibits a downward peak because the aeration/bubble is partly above ℎ line  .This illustrates how the aeration below ℎ line  is added to ℎ line  to yield ℎ impr  .

Fig. 6 .
Fig. 6.Approximation of the liquid level by ℎ bin  for two different thresholds  bin = 0.5 and  bin = 0.8 plotted in the liquid volume fraction field for   = 0.35 m/s (left), 1.87 m/s (middle), and 4.0 m/s (right).

Fig. 7 .
Fig. 7. Visualization of  up and  low for exemplary cases.

Fig. 7
Fig. 7 illustrates the choice of upper and lower thresholds for ℎ line  for four different cases.For each case, a slug structure is marked with a red rectangle, and a plot of the liquid volume fraction of the corresponding time interval is shown on the right.This illustrates the fraction of liquid and gas inside the slug structure.For the case with a superficial gas velocity of 0.62 m/s, four slug structures and one wave at around 21.75 s can be visually identified.If the liquid level of a slug structure fluctuates around  up , e.g., at 22.4 s, then  low prevents double counting of slugs until the liquid level falls below  low .The liquid volume fraction plot on the right hand side shows a slug structure and the wave.It also shows that the fluid phases are well separated.For the 2.5 m/s case, the liquid level does not exceed  up between 21.0 s and 21.5 s.However, the liquid volume fraction plot reveals that the liquid phase is aerated and reaches the top of the pipe.This demonstrates how aeration affects the calculated liquid level and why  up has to be adapted to the superficial gas velocity or slug holdup.However, each slug is unique and may have a different shape and level of aeration.Therefore, a certain  up might be too high to detect all slugs, but too low to also count waves as slugs at the same time.Comparing the liquid volume fraction plots on the right, one observes how aeration increases with higher superficial gas velocities.In the plots of the liquid level time series (left picture), on the other hand,

F
.Webner et al.

Fig. 8 .
Fig.8.Slug frequency   over upper threshold  up for five selected test cases with and without utilization of lower threshold  low .The circles mark the chosen  up for these cases and the respective slug frequencies.

Fig. 9 .Fig. 10 .
Fig.9.Slug frequency   over lower threshold  low for five selected cases.The circles mark the chosen  low for these cases and their respective slug frequencies.

F
.Webner et al.

F
.Webner et al.

Fig. 12 .
Fig. 12.Comparison of evaluation methods with slug frequency prediction methods from literature.

F
.Webner et al.

Fig. B. 1 .
Fig. B.1.Slug frequency predictions from literature.The markers indicate the test cases considered in this paper.
836 + 2.75 sin 0.25 ()) s −1 ,(B.4)where  is the inclination angle.For an inclination of 0 • , as discussed in this paper, Eq. (B.4) yields: (2009) analyzed 230 data points published in literature (where the internal diameters ranged between 2.5 cm and 20.3 cm) to develop a new slug frequency prediction method:

Table 1
Properties of the fluids used in the numerical simulation.
• ] 72 Fig. 1.Simulated test points in the flow pattern map by

Table 2
Spatial discretization schemes used in the numerical simulation.

Table 3
Parameters used for the calculation of pwelch.

Table 4
Respective impact of changing  low and  up by 0.05 on slug frequency.
found the following correlation: developed a new correlation based on the collected data: Eq. (B.7) is for horizontal and turbulent flow.For laminar flow, typically fluids with high viscosity, and/or inclined pipes, correction factors can be applied.Comparison of the slug frequency prediction methods.Fig. B.1 shows the described slug frequency prediction methods from literature for horizontal flow at a superficial liquid velocity of 1.87 m/s and an internal diameter of 0.097 m.The slug frequency   is plotted against superficial gas velocity   .The markers refer to the flow conditions considered in this paper.