Performance enhancement of axial concurrent liquid–liquid hydrocyclone separator through optimization of the swirler vane angle

Vane angle configuration is considerably affecting the internal flow behavior and separation performance of a concurrent axial inlet liquid–liquid hydrocyclone. This study was carried out to improve the design of the swirl generator by optimizing the vane’s deflection angle in an oil/water axial inlet hydrocyclone separator. Angles ranging from 37° to 75° were examined at various operational conditions, including mixture temperature, mixture flow rate, and water-to-oil ratio. Two analysis techniques have been coupled to achieve the aim. First, design of experiment by the response surface method was utilized to generate a combination of run/boundary conditions of swirler vane angles, inlet mixture temperatures, flow rates, and concentrations. The obtained 15 run/boundary conditions were adopted as cases for computational fluid dynamics simulation to determine the separation efficiency, tangential velocity and pressure drop of each case using ANSYS Fluent software. The optimization results show that the swirl generator with a 45° deflection angle generated slightly higher tangential velocity compared with higher and lower vane deflection angles. The separation efficiency obtained by using the 45° swirl generator was higher than other angles, in spite that the turbulence intensity is slightly higher at 45° compared to other vane angles.


Introduction
Numerous oil wells are subjected to the waterflooding technique to maintain reservoir pressure and exploit it to the maximum, especially in shallow offshore and deep water. Even onshore wells are experiencing such a problem (Zhao et al. 2017;Wang et al. 2019). As oil fields become mature, the production of water cut increases considerably under the effect of conning phenomena. Gravity separators are widely used to separate the oil/water, at the surface. On other hand, cyclones have been used extensively to separate solid/liquid and solid/gas cyclones in industries. Liquid/liquid hydrocyclone (LLHC) technology is adopted to tackle the problem for oil/water separation technique. However, due to the small density difference between phases, the application of hydrocyclones to the separation of two immiscible liquids, such as oil and water, is uncommon (Lu et al. 2019). Thus, a highly compact separator in terms of size and performance should be developed to replace the existing gravity-based bulky separation technology.
LLHCs have been used as a downhole oil-water separator to treat the extra produced water. Two types of LLHC separation technology exist, and the fundamental physics behind each separation process is centrifugal force. This force depends on the input and output in the cyclone either tangentially or axially and the outputs from the top or bottom, which is referred to as countercurrent or concurrent axial flow, respectively. Low pressure drops, low droplet breakup, and a stable vortex are several of the advantages of axial concurrent flow cyclones over countercurrent reverse flow cyclones. Thus, axial concurrent flow cyclones have a wide utilization scope in the future of the oil industry both at surface and downhole production intervals (Kitoh 1991;Dohnal and Hajek 2016). In an axial flow separator, centrifugal force is achieved by a liquid mixture flowing through a stationary swirl generator placed inside the separator, and the fluid moves concurrently downward.
Several studies have performed modeling, characterization, and optimization of accessible reverse flow in tangential inlet cyclones, whereas only a few studies have dealt with axial inlet concurrent hydrocyclones. The primary purpose of the swirl generator in axial concurrent flow is to induce high tangential velocity for achieving good centrifugal flow. An excellent axial velocity must be provided at the inlet to achieve decent fluid deflection. A jet effect is observed at the end of the trailing edge of the swirl generator, and high tangential velocity is obtained. The flow should not create any blockage nor flow deformation. However, no completed design criteria are available for the blade profiles in axial LLHC. Dirkzwager (1996) proposed the use of axial concurrent cyclones for oil-water separation. He found that a swirl decays downstream, which hampers separation. For an efficient system, the swirling strength must be sufficiently high but should not cause shear-induced droplet breakups. However, Dirkzwager did not disclose any information on swirl and nozzle designs, and this lack of information limits 3D numerical work for in-depth numerical research. Shi et al. (2012) experimented with a cylindrical-conical cyclone with holes in the conical section at 8% volumetric oil content by using a slit-type outlet for the water stream. Vanes (40° angle) constructed from semicircular plates were determined to be excellent in coalescing the oil droplets. The tested concurrent design had a no-reverse-flow region anywhere within the body as opposed to the standard version. However, the researchers only measured the axial velocity in one location and regarded the unit as "reverse flow free." This claim contradicts the results of other experiments, such as those of Slot et al. (2011), who found an evident reverse flow.
Detailed computational fluid dynamics (CFD) investigation of flow fields is highly useful in studying the swirl generator effect, but this is impossible due to the absence of blade profiles. Stone (2007) investigated a concurrent cyclone using a tangential inlet cylindrical cyclone at various oil volume fractions (low, moderate, and high). The highest efficiency of 80% was achieved at 75% water cut by using a 3× diameter of cylindrical length. Poor performance was observed at high oil loading of 65% and above. These results might have been caused by intense droplet breakups due to a tangential inlet that could have hampered the separation efficiency. Inlet velocity was set to 1 m/s, which is considerably high for a 210 mm long and 60 mm wide cyclone. However, Stone (2007) could have counteracted swirl decay by reducing the dimensions but probably at the cost of droplet breakups, which eventually reduce efficiency. He suggested that a conical body can be used to reduce swirl decay and increase residence time. This suggestion has been considered in the current investigation. The axial LLHC investigated in the present article has conical body.
Swirl generator studies were conducted by Rocha et al. (2009). The researchers studied the effects of utilizing a pointy cone swirl generator tail in the flow using a cylindrical body. The design effectively reduced the reverse flow and increased the tangential velocity in the cyclone probably because of the smooth flow transition. Rocha et al. (2009) suggested that the best cone angle for reducing vortex breakdown is related to the guide vane swirl angle. However, no geometrical data on the swirl generator were provided. Van Campen (2014) used a cylindrical cyclone with 100 mm diameter and 1700 mm total length. Three swirl generators "weak," "strong," and "large" were used, and they differed in terms of flow deflection and hub radius. An annular reverse flow was observed in the experiment, and it was eliminated because it prevented droplets from migrating to the center. The swirl strength also decayed axially due to wall friction that radically deteriorated efficiency. A conical body should have been used instead to allow increased centrifugal forces to act on the droplets, as Van Campen (2014) concluded. He did not disclose enough information on swirl design. CFD has also been used to design the blade of an axial flow hydrocyclone by Zhang et al. (2014). The fluid behavior inside the hydrocyclone can be analyzed with the change in vane angle profiles. "Best vane angle profile can be optimized to achieve higher performance" is what Zhang et al. (2014) recommended.
Modification and optimization of swirl generators in axial inlet concurrent hydrocyclones were discussed by Dirkzwager (1996), Slot et al. (2011), andVerlaan (1991). In their researches, most of the information related to the modeling of the swirl generator was based on some undesirable assumptions that were either not complete or not solidly satisfactory. Varying the vane exit angle affects the exit velocities and swirling efficiency. Thus far, no studies have shown the importance of using the response surface method (RSM) for the optimization of the swirler vane angles.
The primary motivation of the current research is to develop and optimize swirl generators of an axial inlet concurrent LLHC. The analysis procedures used in the current study, namely RSM for optimization and CFD for simulation, are neatly coupled. RMS has been utilized to produce the run sets of design and operation parameters. The run sets, then simulated in CFD environment, where various system performance parameters, such as separation efficiency, pressure drops, and flow structure parameters, such as tangential velocity, turbulence intensity, and static pressure distributions have been predicted at various operational conditions. Additionally, the tangential velocities obtained from the CFD analysis have been incorporated into the design of experiments (DOE) to optimize the vane angles. The selected process variables in the modeling and optimization are four vane angles, five inlet mixture velocities, and five inlet mixture temperatures. The ranges of variables are close to the actual field state of oil production and separation.

Methods and materials
An axial inlet hydrocyclone operates based on the principle of centrifugal force. The flow is generated axially at the top (inlet) and passes through the swirl generator to the bottom. The swirl generator placed inside the hydrocyclone body helps generate the swirling flow, which in turn separates fluid mixtures with a dissimilar density. The vanes deflection angles are vital parameter on the separation performance of the axial inlet concurrent LLHC.
The key procedures implemented in this study are RSM and numerical simulation, which are discussed in detail in "Response surface method" and "Numerical method" sections, respectively. To optimize the experiments, only one parameter is adjusted at a time while keeping the remaining constraints constant. Predicting the interaction between constraints is time-consuming and difficult. To manage this difficulty, RSM is utilized together with numerical simulation to investigate the response and interaction of three operating parameters. The coupling of the RSM and CFD is utilized according to the simple flowchart outlined in Fig. 1.

Response surface method
DOE is widely used in nearly all industries to understand the behavior of parameters. It helps reduce the number of unwanted experimental runs and thus saves time and material resources. Response surface model design is useful in investigating the nature of the relationship between responses and factors, especially when the relationships are nonlinear, in comparison with determining the essential factors. This design helps model the quadratic effect of factors and is well suited for predictive modeling and optimization. RSM is useful in planning, conveying, developing, and analyzing new scientific learning and products. It is also efficient in improving current studies and products. RSM is applied in nearly all areas, such as industrial, biological, chemical, food, and engineering sciences (Oehlert 2000;Bradley 2007).
The first objective of RSM is to determine the optimum response in the presence of multiple responses, identifying the compromise optimum that does not optimize only one response is crucial. The second objective of RSM is to determine how the response varies in a specified direction by adjusting the preferred variables as outlined by Aanchal and Knika (2016). This study used the standard design portfolio of the Design of Experiment Centre's composite design. A surface design is created for investigating the relationship between process variable, such as inlet velocity, temperature, and vane angle response to obtain the exit velocities from the optimum vane angle. By inputting the minimum and maximum values of each variable, data sets/boundary conditions of different combinations are generated. Minimum and maximum values of variables provided to the DOE for the vane angle, inlet mixture velocity, and mixture temperature are provided in Table 1.
Total of 15 combinations are generated by DOE, as shown in Table 2. By default, DOE generated some cases with values less than the minimum and higher than the maximum factors. Each combination suggested by the DOE is representing input boundary conditions to the CFD procedure. Each one of the 15 boundary condition generated by the DOE is inputted to the CFD simulation, and the exit tangential velocities are obtained. The boundary conditions and the resulted tangential velocities are then put back to

Numerical method
CFD is a cost-effective means of understanding complex fluid behavior and its effect on various geometrical and operating parameters. With the aid of powerful computers, the issue of solving strongly coupled, nonlinear partial differential equations of mass and momentum conservation associated with fluid behavior has become easy.

Modeling and meshing
The geometry requirement of the axial inlet concurrent LLHCs is that the total length should be 1330 mm. The geometrical dimensions of each section in the hydrocyclone are shown in Fig. 2. The inlet flow direction from top is a cylindrical shape with diameter, D = 40 mm and length, L 1 = 70 mm. Down the inlet, a stationery swirl generator is placed inside the hydrocyclone to induce swirling flow, as shown in Fig. 2. The upper part of the swirler zone, the nose part, is diverged from D = 40 mm to D 1 = 60 mm with 80 mm length and continue as 60-mm-diameter cylindrical section for L 2 = 610 mm below the nose point. Then, the passage starts to converge from D 1 = 60 mm to D 2 = 40 mm forming the conical part which has a length of L 3 = 430 mm. Then, the outlet section, of length, L 4 = 200 mm, converges D 3 = 20 mm, which is the diameter of light phase outlet, LPO. The swirler has three parts, nose, central ring, and tail, as shown in Fig. 3, with total length of 240 mm. Each part has 80 mm length. The nose is the top part of the swirl generator geometry, which has an elliptical shape, followed by the vane and inner ring. The bottom radius of the nose is like the diameter of the central ring, 42 mm. The cross section at the vane exit is considered as the reference to the elevations, i.e., Z = 0. Above this, the elevation is negative, and below that, the elevation is positive. The nose guides the inlet oil/water mixture toward the annular region of the swirler and thus increases the mixture flow velocity toward. Seven asymmetrical vanes have been produced separately, using CNC machining, and attached to the central cylinder of the central ring. The bottom part of the swirl generator geometry is the tail, and it is elliptical shape. The tip of the tail is not sharp; instead, it is rounded to allow smooth transition of the flow and to reduce any fluid deformation. The model, shown in Fig. 3, is created using CATIA V5 and fabricated as physical model using CNC machine. The mesh and cell discretization is performed in this study by a hyper mesh software. A hybrid meshing scheme, shown in Fig. 4, is selected. Hybrid mesh integrates the structured meshes, used in the regular geometries, and the unstructured meshes, used in the swirler complex flow region, in an efficient manner. Structural hexahedral elements are necessary to reveal the flow distribution beneath the swirl generator. Hence, the meshing is divided into three sections of the LLHC separator of this study. First section is the top cylindrical inlet section meshed with the structural elements. Second section is the unstructured elements used throughout the swirl generator section. Third section is the structural elements below the swirl generator sections. In all three sections, hexahedral elements are adopted six faces sharing with neighbors tends to produce more accurate results of the flow. Mesh quality is gauged by skewness, which is maintained as 0.6 for the meshing. A mesh with total of 980,000 elements has been adopted for the study, and the solution is verified for grid independent.

Multiphase model and boundary conditions
Navier-Stokes equation is solved using the finite volume method to determine the exit tangential velocity from different swirl generators. The density and viscosity of brine at 25 °C were 1067 kg/m 3 and 0.00089 kg/m s, respectively, whereas those of oil were 835 kg/m 3 and 0.01506 kg/m s, respectively. Varying velocity inlet conditions and varying velocities from 0.5 to 2.5 m/s were fed to each DOEgenerated case at different vane angles a combination of runs is generated by the DOE to determine the effect of the inlet flow rate and temperature with varying vane angles. ANSYS Fluent is used to run the boundary conditions for the numerical simulation where, no slip at the wall, time regime is steady state. Multiphase model is applied with SIMPLE pressure-velocity coupling algorithm and second-order upwind discretization scheme. The least square cell-based interpolation method is used. The outlets are imposed with an outflow boundary conditions.
The flow within LLHC is complex and swirls in different directions. Therefore, a turbulent model that can accurately determine the anisotropic turbulence within the cyclone is required. Reynolds' stress turbulence model is a secondorder turbulent closure that can capture anisotropy, has a high streamline curvature, can change the strain rate quickly, and allows for severe swirling conditions to achieve realistic simulations of cyclone flows (Utikar et al. 2010). Thus, it has been chosen for this study. The model can be expressed, after Liu et al. (2010), Cui et al. (2014), as follows: where P ij is the stress production term defined in Eq. 2, DT ij is the turbulent diffusion term defined in Eq. 3, Φ ij is the pressure-strain term, ε ij is the viscosity diffusion term, and F ij is the rotation production term. They are defined as follows:

Validation of the CFD procedure
To ensure an accurate numerical simulation, the simulation results are usually associated with those of similar previous validated numerical works. The simulation results in this study were compared with Slot et al. (2011) numerical results. The tangential velocity at 0.44 m distance from the vane angle toward the bottom is depicted in Fig. 5. The trend of predicted tangential velocity profiles in both works, Slot et al. (2011) and the current work, is identical with mean deviation below 25%, which denotes an acceptable agreement. The minor difference is attributed to the difference in the meshing procedure. (1)

Fig. 3: 45° vane angle created with respect to the horizontal axis (left) and modeled swirler (right)
Furthermore, comparison of experimental results and numerical results of 80:20 water-to-oil ratio at 80 °C and 90:10 is shown in Figs. 6 and 7, respectively. The results include the separation efficiency versus the flow rate obtained from both experimental measurements and numerical simulation results. The results of numerical simulations showed good agreement with experiment. However, an over prediction of nearly 9% of separation efficiency from the experimental results is observed. The validation confirms that the simulation procedure adopted in the current research is able to predict the system performance and the flow behavior with good accuracy.

Results and discussion
The flow distribution inside the LLHC has been determined by examining the separation efficiency, pressure drop, radial distributions of the tangential velocity profiles, turbulence intensity and static pressure. Qualitative analysis, by velocity contour plots, and quantitative analysis, by variables prediction and plotting, have been utilized to analyze the results. Optimization by RSM is also presented and discussed. Location of Z = 90 mm has been identified to compare the tangential velocity, turbulence intensity, and static pressure in radial distributions for 45° and 68° vanes angles.

Analysis of separation efficiency
The flow rate vs separation efficiency measured for the best and least deflection vane angles is shown in Fig. 8. The flow rate vs separation efficiency is high for 45° vane angel compared to the 68° vane angle. It is examined that the separation of oil-water mixture requires higher tangential velocity inside the LLHC to have a good efficiency. For a lesser flow rate, the efficiency is minimum. Larger flow rates of the mixture start hike in the separation efficiency. It is postulated to be the lesser tangential velocity makes the 68° vane angle for a reduced separation efficiency. However, with the increase in flow rate, even in the case of 68° vane angle, the efficiency is increasing considerably.

Analysis of pressure drop
The pressure drop has been studied for all the 15 cases generated by DOE, listed in Table 2, and simulated by CFD. The pressure drop results are shown in Fig. 9. A trend of increasing the pressure drop from high vane angle toward low vane angle. That is postulated to high tangential velocity at the exit of the vanes causing larger shear losses with the cyclone wall and the tail body of the swirler. The larger shear i.e., higher oil droplets breaking, is an added agitation for larger instability of the mixture and lesser oil/water separation. The pressure drop is acute for the vane angle lesser than 45°, indicates high instability of flow field inside the LLHC, and leads to the consumption of more energy.

Analysis of tangential velocity
The nature of the flow in a cyclone separator is mainly influenced by its tangential velocity component. The velocity contours presented in the longitudinal cross section, shown in Fig. 10, indicate that the maximum velocity is generated at exit from the vanes. However, the tangential velocity developed to the maximum beneath the swirl generator, and a reduction in tangential velocity is observed when the flow progressed toward the bottom. This trend of the result is in good agreement with the results of Oakman (2015) and Noroozi and Hashemabadi (2011).
To quantify the difference in the produced tangential velocity, Fig. 11 is produced. The figure displays the tangential velocity distributions results across the cone diameter, at elevation Z = 90 mm, produced by a 45° and 68° vane angles. The maximum predicted tangential velocity, of 7.06 m/s, is thrice the average inlet velocity for a 45° swirl generator. Meanwhile, less maximum tangential velocity, of 4.17 m/s, is noted for the 68° deflection angle. Two zones of vortex flow existed inside the LLHC can be realized. Near the wall, low tangential velocity was initially present. The trend of increasing and reaching the maximum velocity depicts a potential vortex flow at a location slightly away from the wall. Then, the magnitude of the tangential velocity starts to reduce gradually toward the central zone. At the longitudinal central line, i.e., at r = 0, the tangential velocity is almost zero, which is similar to the findings of Bergstrom and Vomhoff (2007). When moving toward the center, the tangential velocity distributions become minimum, which resembled solid body rotations, and their widths are proportional to the swirl intensity distribution along the axis of the cyclone. A reduction in the tangential velocity distributions toward the axial distances is observed to be faster for the case of 68° vane angle. In the solid body rotation center, the slope of the tangential velocity distributions is maximized and moved toward the core region; then, it increased toward the side of the wall (Al-Kayiem et al. 2014;Mokni et al. 2015).

Analysis of turbulence intensity
Depending on the magnitude of inlet feed velocity, turbulent fluctuations become unpredictable and complicated most of the time. This condition is postulated to be due to the mixing and colliding of streams, thus causing recirculation and flow reversal. The turbulence inside a hydrocyclone is proportional to the difference in velocity gradients as realized in   Fig. 7 Comparison between experimental measurements and CFD prediction of separation efficiency results of, 90:10 water-to-oil ratio at 80 °C, for different flow rates Gomez (2002) experiments. The turbulence intensity, T.I., also often named as turbulence level, is defined as the ratio of the root-mean-square of the velocity fluctuations, u ′ to the mean flow velocity,U and can be predicted using Eqs. (4) to (6), (see Russo and Basse 2016;Wang et al. 2007).
The root-mean-square of the velocity fluctuation is the mean of 3D components, x, y, and z.

And the mean velocity, as
The turbulence intensities generated from vane angles 45° and 68° vanes were plotted to examine the distribution of turbulence across the radial distance at Z = 90 mm and are shown in Fig. 12. The turbulence intensities increased considerably beneath the swirl generator, which could be due to the rapid flow development immediately after the hemispherical bottom section of the swirl generator (tail). In Fig. 12, the line with circular symbols denotes the turbulence intensity distributions for the 45° vanes, and the line with triangular symbols denotes the turbulence intensity distributions for the 68° vanes. Turbulence intensity is at the maximum for the 45° vane angle, and the trend of maximum turbulence intensity toward the wall decreased in the direction of the center in both vane angles, 45° and 68°. The high turbulence at the wall is believed to be due to the flow being uneven, which may be because of the high tangential velocity. However, low turbulence intensity is noted from the 68° vane angle beneath the swirl generator and is shown in red color lines.   Figure 13 shows a comparison of the static pressure distribution across the radial distance of Z = 90 mm. In Fig. 13, the line with circular symbols indicates the static pressure depicted from the 45° vane angle. The trend of increased pressure near the wall and toward the core is reduced considerably. This result can be correlated with Fig. 9, in which high tangential velocity is observed at the same point. Hence, the high energy consumption of high static pressure is due to high tangential flow near the wall. This sudden pressure drop supports the presence of backflow due to the formation of a high swirling flow or high energy consumption just beneath the swirl generator. The static pressure distribution from the 68° vane angle is at the minimum and showed minimal energy consumption due to the reduced swirling flow. Furthermore, the minimum tangential velocity for the 68° vane angle at the same location justifies the above statement.

Optimization analysis of response surface method
The obtained tangential velocities from the numerical simulation as per DOE cases are shown in Table 3. The optimum responses for vane angle, especially for the combinations of inlet velocity and temperature, have been obtained with the aid of the response optimizer in DOE software.
The maximum tangential velocity, 7.32 m/s, imposed by the vanes is exhibited by the 45° vane angle with an inlet velocity of 2.0 m/s. The least resulting tangential velocity is noted with 0.5 m/s inlet velocity and 68° vane angle. At 2.0 m/s inlet velocity, the exit tangential velocities in 56° and 68° vane angles are comparatively lower than the case of 45° vane angles. The change in the response to a direction when altering a specific variable could be captured by scheming surface and contour plots, as presented in Fig. 14 (top). The surface and contour plots are supporting the decision on the optimization. These plots show only the response of two factors, which are inlet velocity and vane angle, at a time by keeping the third factor (temperature) constant.
An interaction plot, shown in the RSM analysis in Fig. 14  (bottom), is studied to justify the optimized design of the vanes deflection angle. The plot is used to examine the difference between a "level means" of more than one factor, simultaneously. An interaction effect occurs when diverse levels of a factor influence the response differently. An interaction plot displays means for the levels of one factor on the x axis and a distinct line for each level of another factor. The lines, whether parallel (no interaction) or nonparallel (interaction), are analyzed to determine the nature of the relationship and its influence on the relationship between the factors and response (Granato and Calado 2013). Interaction plots, shown in Fig. 14, indicating that 45° vanes angle and 2.0 m/s are resulting in the highest tangential velocity, value of 7.32 m/s. In addition, ANOVA analysis is performed, to determine the significant parameters. Commonly, the significant value is determined by the P value, which must be less than 0.05 (Panagiotis and Nikolas 2015). ANOVA analysis shows that the inlet velocity, inlet velocity square, and the interaction of inlet velocity with vane angle are significant. Thus, it has been concluded that inlet velocity of the mixture and vane deflection angles are the most significant parameters in optimizing a swirl generator.
The regression analysis above led to the establishment of a mathematical equation for the statistical relation of inlet velocity and vane angle to the response variable, which is exit velocity.
where Y is the response variable (tangential velocity), X 1 is the inlet velocity, and X 2 is the vane angle. The temperature was not included as its P value is very high which indicates that it is not significant influencing variable.

Conclusion
An effective and successful analysis and optimization procedure are developed in the current study by integration of DOE, CFD, and RSM techniques. The developed procedure has been utilized to reveal the optimum vane deflection angle in axial concurrent flow oil/water separator. Optimization results by RSM showed that the swirl generator with a 45° deflection angle and an inlet mixture velocity of 2 m/s generated high tangential velocity, whereas that with a 68° deflection angle and 2 m/s inlet mixture velocity resulted in lower tangential velocity. Aside from tangential velocity, turbulence intensity was also investigated. The turbulent intensity generated by the 45° swirl generator was higher than that produced by the 68° swirl generator. The design of swirler with 45° vane deflection angle has resulted in highest separation efficiency of the axial inlet concurrent liquid-liquid hydrocyclone separator, and it is the optimum. Hence, 45° is optimum vane deflection angle to design the swirler of axial flow liquid-liquid hydrocyclone separator. A mathematical model of a swirl generator for an axial inlet hydrocyclone was developed and statistically established based on RSM.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creat iveco mmons .org/licen ses/by/4.0/.