Numerical investigation of the effect of viscosity in a multistage electric submersible pump

ABSTRACT Electric submersible pump (ESP) systems are commonly used as an artificial lift technique by the petroleum industry. Operations of ESPs in oil wells are subjected to performance degradation due to the effect of oil viscosity. To understand this effect a numerical study to simulate the flow in three stages of a multistage mixed-flow type ESP operating with a wide range of fluid viscosities, flow rates, and rotational speeds was conducted. The problem was solved by using a commercial computational fluid dynamics (CFD) software. The numerical model was validated with experimental head curves from the literature at different viscosities and rotational speeds available for the same ESP model used in this study, and good agreement was found. Performance degradation was evaluated by analyzing the effect of viscosity on head and flow rate. In addition, a flow field analysis to compare the flow behavior when the pump operates at different viscosities was carried out. The interaction between stages was also analyzed, and the influence of a previous stage on the upstream flow was evidenced. The flow field was analyzed at a curved surface that follows the complex mixed-flow geometry of the stages. CFD proved to be useful for exploring this kind of feature, a task whose accomplishment by means of experimental methods is not trivial. Such analysis helps to understand the flow pattern behind head and flow rate degradation when the Reynolds number is decreased. The results from this work are helpful as they provide a basis to estimate performance degradation for general scenarios.


Introduction
Artificial lift techniques are used to increase the production of a given oil well. Electric submersible pump (ESP) systems are a common example of that technique and are used in both onshore and offshore installations. The main component of an ESP system is the multistage centrifugal pump, which provides the energy required to lift the fluid from the wellbore to the production facilities, be those latter onshore or offshore.
One of the problems of using ESPs to pump crude oils is related to the performance degradation with respect to the regular operation with water that the pump experiments due to the effect of the fluid viscosity. Once friction losses can increase significantly depending on the viscosity, its influence on performance degradation is twofold since a higher power input is required by the ESP, whilst the pump head and flow rate decrease. Eventually, this combined effect severely decreases the pump hydraulic efficiency.
Performance degradation of centrifugal pumps operating with highly viscous liquids is a subject that has been investigated for years. There are engineering standards and procedures to predict performance degradation, such as the ones proposed by Stepanoff (1967), Gülich (2010), andHydraulic Institute (2015). However, in almost all cases the range of specific speeds for which such methods are valid encompasses only radial impeller pumps, whereas the majority of ESPs used in offshore wells have mixed-type impellers.
In addition, there is a lack of information about the flow behavior of highly viscous liquids inside centrifugal pumps, especially for multistage mixed-flow types. Amaral (2007) and Amaral, Estevam, and França (2009) studied the effect of viscosity on the flow in different centrifugal pumps, including a mixed-flow type ESP. The authors experimentally evaluated performance parameters such as head and efficiency for a wide range of viscosities and developed an algebraic model to predict the performance of each pump in those situations. However, no thorough attention was given to the performance degradation itself or to associated flow field phenomena. Maitelli, Bezerra, and da Mata (2010) numerically studied the flow inside one stage of a mixed-flow type ESP. The authors used a non-structured numerical grid to simulate one stage of the ESP. Initially, the impeller was analyzed separately, and the resulting head values were found to be higher than the values taken from the manufacturer's data sheet. When the diffuser was added downstream of the impeller, the head curve showed better agreement with the curve provided by the manufacturer. Vieira, Siqueira, Bueno, Morales, and Estevam (2015) developed an analytical model to predict head curves of an ESP operating with oils. The author analyzed a combination of several analytical models to predict the pressure loss in different parts of the centrifugal pump. Sirino, Stel, and Morales (2013) numerically studied the performance of a single stage of an ESP. Large instabilities in the water head curves were observed, which were explained by the authors as a consequence of using only a single stage in the simulations. Stel et al. (2014) numerically investigated the performance of three stages of an ESP. Some flow aspects between simulations with water and a more viscous fluid were compared. The authors observed that the head curve instabilities found by Sirino et al. (2013) were greatly reduced with the three-stage numerical model. The performance throughout each stage was compared, showing that there are significant differences in the pressure gain provided by each stage. In addition, the authors evaluated the performance degradation for a constant rotational Reynolds number. The authors showed that for different conditions but equal rotational Reynolds number, one can observe similarities in the performance curves if proper dimensionless parameters are compared. Recently, Stel, Sirino, Ponce, Chiva, and Morales (2015) developed a numerical study on the flow inside multiple stages of an ESP operating with water at several flow rate conditions. In this study, the authors verified that the number of stages considered in the simulations has an influence on the average performance results. In general, the first stage shows the greatest difference with respect to the pressure rise of the following stages, whereas from the third stage onwards the pressure rise value of each stage changes only slightly. In addition, the authors explored the flow field inside the ESP at different flow rates. They observed reverse flow zones for part-load conditions in the diffuser, which restricted the stage net outflow area. Different flow patterns were verified inside the ESP depending on the flow rate operation.
Yet, the flow behavior inside the ESP operating with oil is not clearly understood. Transition from turbulent to laminar flow regime may occur when the ESP is operating with oil. In addition, the flow field from diffuser to the next impeller was not explored for fluids with different viscosities. Analyzing flow field features in a mixed-flow multistage pump stage by experimental methods is quite difficult, due to the curved geometry of the hydraulic channels.
Numerical tools, on the other hand, can be very helpful to handle complex geometries, at the same time providing flexibility when testing a wide range of operating conditions. Review of the literature shows a vast number of works dealing with pump flow analysis using CFD. For instance, Shi and Tsukamoto (2001) used CFD to investigate pressure fluctuations caused by the interaction between the impeller and a vanned diffuser. Yedidah (2008) presented a critical analysis on the use of CFD when studying flows in pumps. Feng, Benra, and Dohmen (2010) used a commercial CFD software to test the capability of several turbulence models in predicting both local flow field features and performance data. Feng, Benra, and Dohmen (2007), Stel et al. (2013), Liu, Wang, Wang, Huang, and Jiang (2014), Tan, He, Liu, Dong, and Wu (2016) and Yao, Yang, and Wang (2016) also employed CFD to study several flow aspects in pumps. However, most numerical works on the subject assume water as the working fluid, being the effect of viscosity rarely taken into account.
Based on the aforementioned considerations, this article presents a numerical study of the flow in three stages of a mixed-type ESP. The influence of viscosity on the performance and flow field is investigated in this study. In addition, some aspects of the ESP operating with fluids of different viscosities are clarified, such as the interaction between the diffuser and the impeller. Understanding all of these features is important to comprehend the phenomena associated with performance degradation and, in extension, to provide better estimations for the performance of pumps under those conditions.

Pump geometry
The Schlumberger Reda TM GN7000 540 series is the multistage ESP model used in this numerical study. Amaral (2007) and Amaral et al. (2009) used a three-stage mount of the same ESP model operating with fluids of different viscosities to measure performance parameters. Stel et al. (2015) also used a three-stage configuration of the same ESP in their numerical study, which is assumed again in this work. The impeller and the diffuser geometries were digitized with the aid of 3D scanning of a real Schlumberger Reda TM GN7000 540 ESP stage and later reproduced numerically in ANSYS TurboGrid TM .
The aforementioned ESP has a mixed-flow impeller with seven blades and a diffuser with seven vanes.   Figure 1 presents a picture of a single stage of an ESP, with indications of the main diameters and positions of the hub, shroud, blade, vane, inlet, and outlet of the stage. A virtual midspan surface, which will be later used for some analyses, is also schematically represented. This surface follows the mixed-type geometry all the way from inlet to outlet, being located halfway between the hub and the shroud. Table 1 shows the dimensions presented in Figure 1, where the subscripts 1, 2, 3 and 4 refer to the impeller inlet, impeller outlet, diffuser inlet, and diffuser outlet, respectively. Subscripts i and o correspond to the inner and outer diameters depicted in Figure 1. Other information about the ESP dimensions such as the blade angles and thicknesses can be found in Stel et al. (2015).
According to the manufacturer's data sheet, when this ESP is operating with water at a rotational speed of n des = 3500 rpm, it should deliver Q des = 1.360·10 −2 m 3 /s and H des = 9.6 m per stage, at its best efficiency point (BEP). This gives a specific speed of ω s = des · √ Q des /(H des · g) 0.75 = 1.41 (with rad/s, m 3 /s, m and m/s 2 ), or n q = n des · √ Q des /(H des ) 0.75 = 74.6 (with rpm, m 3 /s and m).

Numerical model
Turbulent, transient, single-phase, isothermal, and incompressible flow in the ESP is modeled through the unsteady Reynolds-averaged Navier-Stokes equations, usually referred as U-RANS equations. They can be represented in a general form as: where, ρ is density, t is time, and − → V is the Reynoldsaveraged velocity vector. In the continuity equation, φ = 1 and = S = 0. The momentum equations assumes φ = − → V , = μ + μ t and S = −∇p, where ∇p is the pressure gradient, μ is the fluid dynamic viscosity and μ t is the turbulent or eddy viscosity, a term that arises on the grounds of the Boussinesq hypothesis and has to be modeled. For a rotating frame of reference, an additional source term accounts for the Coriolis and the centrifugal effects: where the first and second terms are due to the Coriolis and centrifugal effects, respectively. The shear stress transport (SST) turbulence model (Menter, 1994) has been chosen in this study. This model is not based on wall functions and it is said to handle separation flow generally better than other two-equation turbulence models (ANSYS, 2012). In any case, Stel et al. (2015) tested different turbulence models to evaluate performance parameters for simulations in a single stage of the same ESP studied in this work and found no significant differences in favor of any of them. In addition, the SST turbulence model formulation can be coupled with the transition gamma-theta model (Langtry & Menter, 2005), which is capable of estimating laminar to turbulent flow transition. This is an interesting feature because the operation of the ESP with oils may be subject to different flow regimes, as will be shown later. The gamma-theta model is based on two transport equations: one for the intermittency, γ , and one for the momentum thickness Reynolds number, Re θ . The former variable is used to trigger transition, whilst the latter is solved to set the transition onset criteria (ANSYS, 2012). This model has been extensively validated with experiments for a range of transitional flows (including turbomachinery applications) in works such as Menter et al. (2004), Langtry et al. (2004), Langtry and Menter (2005), Völker (2006), andLangtry (2006).

Numerical implementation
The U-RANS equations are numerically solved by using the commercial fluid dynamics software ANSYS CFX Release 14.5. For turbomachinery problems, it uses a multi-block technique to solve each part of the whole domain as a separated sub-domain. Several numerical aspects of this procedure were presented in detail by Stel et al. (2015). Figure 2 shows the numerical mesh at the hub, shroud and the blades/vanes of one stage and part of the intake pipe. The mesh was built entirely of hexahedral elements using both the ANSYS ICEM CFD TM and ANSYS TurboGrid TM programs. Mesh distribution, especially close to the pump walls, was set in accordance with the requirements of the SST turbulence model, as explained in detail by Stel et al. (2015). A mesh sensitivity test was conducted by Stel et al. (2015) for the ESP operating with water at a part-load flow rate, a case considered very unstable due to the flow instabilities presented in this condition. Further details on mesh arrangement and quality parameters can be found in their work. The chosen mesh for a three-stage ESP has close to 3.6 million nodes, where each impeller and diffuser has approximately 630,000 and 558,000 nodes, respectively. In this work, the same mesh was used for simulations with water as the working fluid, with slight adjustments of the wall node distances for simulations with lower Reynolds numbers (i.e. more viscous fluids). Varying the number of grid nodes for simulations with fluids other than water causes negligible effect in the results. Figure 3 shows the numerical domain of the threestage ESP, as well as the main boundary conditions used to solve the problem. A reference total pressure of 0 [Pa] (gauge) is assumed at the inlet. All walls are smooth and no-slip and a specified flow rate is imposed at the outlet. Two pipes are used to push away the inlet and outlet boundary conditions from the first impeller inlet and last diffuser outlet. The pump inlet and outlet are considered to be at the same hydrostatic reference level. Balancing holes, clearances between impeller and diffuser or any gaps between the stages and the ESP housing were all neglected in this study.
The time-step used corresponds to 12 steps per blade passage, or 84 steps per revolution. The actual time step depends on the rotational speed; for example, 2 · 10 −4 s for 3500 rpm. Stel et al. (2015) showed that increasing the number of time steps per blade passage to as much as 20 causes little influence on the pressure rise through the three-stage ESP when considering a simulation at a part-load condition with water as the working fluid. Simulations with fluids that are more viscous than water proved to be less sensitive to this parameter, thus assuming 12 steps per blade passage for every case is sufficient for the purposes of the present work. Other relevant details of the numerical model can also be found in Stel et al. (2015), such as the total number of time steps simulated, how results are time-averaged to avoid influence of the initial condition on the results and the effect of the number of ESP stages assumed on the pressure rise through the pump.

Flow regime inside the ESP
Determining the flow regime inside a centrifugal pump working with oil is quite difficult. There is no Reynolds number definition that one could use as a way to have an indication of the flow regime. One might use, for example, the rotational Reynolds number, Re = ( R 2 2 /ν), or the Reynolds number based on the hydraulic diameter of the channel, Re D h = (VD h /ν). However, there are no established transition limits based on those definitions that could be used to evaluate or classify the flow regime.
At first, the flow regime was evaluated using the Reynolds number based on the hydraulic diameter of the channel. This was considered the most conservative way to evaluate the flow regime. For fluids with viscosities above 214 mm 2 /s, the Re D h value was found to be below 1000 and then laminar flow was expected. However, when the pump operates with a fluid with a 48 mm 2 /s kinematic viscosity, the Re D h values can be both below and above 2000 depending on the flow rate and rotational speed. As this threshold is uncertain even for flows in circular pipes, one must be aware whether the use of turbulence models around this limit is adequate or not. Then, the cases with a 48 mm 2 /s fluid were simulated using the gamma-theta transition model (Langtry & Menter, 2005)   coupled with the SST formulation and implemented in the ANSYS CFX Release 14.5 software. Figure 4 shows the intermittency variable, γ , plotted at the midspan surface of the ESP operating with a 48 mm 2 /s fluid at three different flow rates (at 1800 rpm) and a 214 mm 2 /s fluid at one flow rate and a 3500 rpm rotational speed. For all simulations, 5% turbulence intensity is imposed at the inlet. The γ value ranges between 0 and 1. A 0 value indicates a fully laminar flow, whereas 1 represents a fully turbulent flow.
From Figure 4 one observes that for an operation at a flow rate of 0.13·Q des the flow becomes turbulent at the impeller, and turbulence is transferred to the diffuser. Large separation zones and unstable flow that happens at part-load operation explain the generation of turbulence for such a low flow rate. However, at the discharge pipe, the turbulence is significantly dampened. For the BEP condition (which is equal to 0.44·Q des when operating with fluids with 48 mm 2 /s at 1800 rpm) the flow is less turbulent than the first case due to a great reduction of the unstable separation flow. For 0.64·Q des , which represents an over-load condition at 1800 rpm for 48 mm 2 /s, the flow enters the impeller with significant turbulent intensity, which in turn is dampened throughout the stages. The last figure shows an operation with a 214 mm 2 /s fluid viscosity at 3500 rpm, which is estimated as predominantly laminar by the gamma-theta model, as first expected when calculating the Reynolds number based on the hydraulic diameter. Although not shown here, all cases using water as the operating fluid are predominantly fully turbulent.
For validation purposes, the liquid properties tested by Amaral (2007) and Amaral et al. (2009) were also assumed in the present work. They used several water and glycerol mixtures at various temperatures to investigate the effect of viscosity on the performance of the same ESP studied in this work. The range of liquid viscosity assumed for fluids other than water is consistent with ranges found in actual petroleum artificial lift scenarios. The liquid properties are described in Table 2, together with the expected flow regime in the ESP. For viscosities above 214 mm 2 /s the flow is considered laminar for all rotational speeds and flow rates simulated in this work. For water, the flow is considered as fully turbulent and the SST turbulence model was used. For the 48 mm 2 /s fluid, the flow is turbulent in some cases and laminar in others,

Results
Numerical simulations were conducted for a wide range of fluid viscosities (1-833 mm 2 /s) and rotational speeds (1800-3500 rpm). Several simulations were carried out assuming various flow rates for a given fluid viscosity and rotational speed, in order to obtain curves of head versus flow rate. The pressure head definition (Srinivasan, 2008) is used in this work. Neglecting any elevation between inlet and outlet, the pressure head is calculated as H = (p in − p out )/(ρg), where p in and p out are areaaveraged static pressure values at the first impeller inlet and the third diffuser outlet, respectively. This definition was used primarily to be consistent with Amaral (2007) and Amaral et al. (2009), whose experimental pressure head data were used to validate the numerical results. In addition, all numerical results are time-averaged values obtained from the transient solution of several impeller revolutions (depending on the variable being analyzed), as explained in Stel et al. (2015).

Comparison with experiments
Numerical pressure head has been compared with experimental data from Amaral (2007) and Amaral et al. (2009). The authors performed experiments with a three-stage Schlumberger Reda TM GN7000 ESP for a wide range of viscosities and rotational speeds. Those references describe the experimental setup in detail. Figure 5 shows the comparison of numerical and experimental curves of pressure head per stage as function of normalized flow rate for different viscosities and rotational speeds. In general, numerical results agree well with the experimental data for both water and fluids with higher viscosities. The coefficient of determination, R 2 , was used to evaluate the agreement between numerical and experimental pressure head values, and a value of 0.9707 was found. Discrepancies may in some cases be attributed to a number of factors, like the impossibility of exactly reproducing the pump design numerically, inherent errors associated to the numerical modeling and experimental uncertainties. Hence, the authors understand that the agreement found is acceptable for the purposes of this work. Figure 6 shows curves of pressure head versus normalized flow rate at a 3500 rpm rotational speed for operation with several fluid viscosities. As it can be observed, the performance of the mixed-flow pump deteriorates continuously with viscosity, which is caused by an increase of the friction losses in the hydraulic channels. The restriction of the operational window is an important result that follows from the degradation of the pressure head and flow rate. For instance, the operational range of the ESP decreases by up to Q/Q des = 0.9 when operating with the highest fluid viscosity assumed (833 mm 2 /s). In addition, the best efficiency point is shifted to lower flow rates when the viscosity increases. Below Q/Q des = 0.4, however, the water pressure head curve falls below some viscous operation curves, a result that is probably related to errors from the numerical model since it is not observed in the experimental results of Amaral (2007). Nevertheless, it should be noticed that this trend is still possible when operation with water and fluids with moderate viscosities are compared, as reported in previous works such as Stepanoff (1967) and Li (2000).

Effect of viscosity on pressure head and flow rate
In addition, a saddle-type instability behavior of the head curve (Gülich, 2010) is observed in Figure 6 for operations with viscosities up to 214 mm 2 /s. This instability causes a decrease or a stabilization of the head with a reduction in the flow rate for a certain range, instead of a continuous rise of the head from high to low flow rates. This behavior is observed below BEP for 1, 48 and 214 mm 2 /s fluid viscosities. For higher viscosities, this instability is greatly reduced.
The phenomenon above is related to the separation flow pattern inside the pump channels at part-load operation, and is especially prominent for diffuser pumps with n q > 30 (Gülich, 2010). Figure 7 shows threedimensional streamlines for three fluid viscosities at partload conditions where an instability in the pressure head curve was observed in Figure 6. For operation with water (1 mm 2 /s), for which the instability in the head curve is the most prominent, one can observe two large recirculation zones both close to the impeller suction side and next to the diffuser outlet. The point at which this instability is observed here for operation with water, Q/Q des = 0.75, is in accordance with the range suggested by Gülich (2010) for this type of pump (0.6 < Q/Q des < 0.9). Inspection of the 3D streamlines also shows that whilst the recirculation zone inside the impeller is mainly blade-to-blade oriented, the swirling flow close to the diffuser outlet is enclosed from hub to shroud. This is also in accordance with the flow pattern related to the saddle-type instability discussed by Gülich (2010). For the lower Reynolds number cases compared in Figure 7, the flow inside the diffuser is better oriented with the geometry of the pump, which may explain why the saddle-type behavior for such situations is reduced when viscosity increases. Figure 8 shows curves of normalized pressure head as a function of the rotational Reynolds number for many normalized flow rates. It is possible to observe that the effect of Re on performance is quite low for Q/Q des = 0.15, even at very low Reynolds numbers (that is, highly viscous fluids). However, the effect of viscosity increases significantly with the flow rate. For example, when operating at a Q/Q des = 1.0 a null pressure head is reached when Re is approximately 1000, whereas for  Q/Q des = 1.25 the ESP loses its pumping capability for a rotational Reynolds number higher than 3000. Figure 9 shows the degradation of the best efficiency point with the rotational Reynolds number obtained for a rotational speed of 3500 rpm. The flow rate associated with the best efficiency point for operation with fluids of any viscosity, Q v,bep , has been normalized by the BEP flow rate for water obtained at 3500 rpm, Q w,bep . It can be noticed that the flow rate degradation with respect to operation with water is relatively low as long as the Reynolds number remains high. However, the Q v,bep /Q w,bep ratio decreases steeply for Re < 10 4 , and the BEP flow rate ratio can become as low as 0.4 when Re is approximately 800. Figure 10 illustrates the friction losses as a function of the normalized wall shear stress for different viscosities and rotation speeds at the hubs and blade/vane walls of the second stage of the ESP. The associated rotational Reynolds number of each case is also shown. All flow rates correspond to the best efficiency point in each case. One can observe a continuous increase in the normalized wall shear stress levels as the rotational Reynolds number decreases. In addition, the blade/vane leading and trailing edges are in general the regions where the highest normalized wall shear stress values are found for each case. When an area average of the normalized wall shear stress over all solid surfaces of the three stages of the pump is taken, the results for a), b), c), and d) are 3.10·10 −2 , 13.11·10 −2 , 21.55·10 −2. and 50.75·10 −2 , respectively. This results in an increase of almost 17 times from an operation with water at 3500 rpm to another with an 810 mm 2 /s fluid at 1800 rpm.

Effect of viscosity per stage
In order to evaluate the influence of viscosity in the performance of each stage of the ESP, Figure 11 presents pressure head curves as a function of the normalized flow  rate for four different viscosities at each stage of the ESP operating at a 3500 rpm rotational speed.
As it can be seen, all stages have a distinct pressure head curve behavior, which is influenced by the flow pattern throughout each stage. Particularly, the pressure head curve of the first stage shows a large deviation from the curves of the following stages for operation with each fluid viscosity, which was also observed by Stel et al. (2015) when the ESP operates with water. In addition, instabilities in the head curve, especially in the case of water for the first stage, resemble the saddle-type behavior previously described. However, an increase in viscosity greatly smooths this effect, which is small for the 571 mm 2 /s first-stage curve. The second and third stages have similar pressure head curves for all viscosities, even though the last stage presents a slightly higher pressure head value near the shutoff point for all cases. In general, the first stage shows a higher head with respect to the following stages through great part of the flow rate range.
At strong part-load operation, however, the first-stage pressure head tends to be lower than for the following stages. Nevertheless, such trend is reduced as the viscosity increases, and is not observed for operation with the 571 mm 2 /s fluid. Figure 12 shows vector plots to illustrate the flow field inside the three stages of the ESP for different rotational Reynolds numbers. All cases depicted are calculated for flow rates at the best efficiency point in each case, and the plots are drawn at the midspan surface between hub and shroud of impellers and diffusers. The velocity component resulting from the projection of the velocity vector over the midspan surface is herein referred to as 'Velocity blade-to-blade' (V bb ), as explained by ANSYS (2012) and Stel et al. (2015). Also, V bb is relative to a frame of reference fixed in each pump part, that is, a non-inertial (rotating) frame of reference at the impellers and an inertial (static) one at the diffusers. The blade-to-blade velocity values are time-averaged from transient results using as many impeller revolutions as needed to achieve a statistically independent value. Also, the blade-to-blade velocity is normalized by the impeller tip speed, U 2 = · R 2 .

Effect of viscosity on the flow field
From Figure 12, it can be observed that an increase in viscosity causes significant changes in the flow field inside the ESP stages. Separation zones are found near the second and third impeller trailing edges operating at Re = 10,229 and Re = 3346. Since friction losses increase considerably when operating with high viscosity fluids and these losses are proportional to the flow rate, the best efficiency point is displaced to a lower flow rate with respect to the pump original design value when the Reynolds number decreases. Once separation flow is consistent with part-load operation, one can expect it to Re Ω 716,046 Re Ω 10,229 Re Ω 3346 Re Ω 644 Figure 13. Normalized blade-to-blade velocity contours at a meridional plane through the three stages of the ESP for different rotational Reynolds numbers. occur in those cases. Recirculation zones are also found near the diffuser trailing edge for operations from Re = 716, 046 down to Re = 3346. This may be a result of a too quick transition to a 90°exit angle of the diffuser vane. Also, it can be noticed that the recirculation zones inside the second and third impellers slightly increase from Re = 10, 229 to Re = 3346. However, for operation at the lowest rotational Reynolds number tested (Re = 644) the flow is well oriented everywhere inside the ESP except in the first impeller. Despite the very low flow rate associated with the best efficiency point for this case (0.15·Q des ), the viscous effect is intense enough to smooth separation zones out. This is also observed in the diffuser where recirculation spots were formed in less viscous cases.
The disoriented flow coming from the diffusers is one of the factors responsible for the appearance of large separation zones in the second and third impellers. Figure 13 presents the normalized blade-to-blade velocity contours at a meridional plane through the three stages obtained at the same operating conditions from Figure 12. One can observe that the velocity distributions at the entrance of the second and third stages are different from that at the entrance of the first stage. This is explained by the fact that the flow inside the second and third impellers is influenced by the disoriented flow from a previous diffuser, whilst the first impeller receives a well-oriented flow from the inlet suction pipe. In addition, the velocity magnitude at the impeller in the two final stages is notably higher closer to the hub than to the shroud, a trend that is different from that what is observed in the first stage. This is caused by the curvature of the impeller and may contribute to the formation of the separation zones found close to the trailing edges of the second and third impellers for Re = 10,229 and Re = 3346. Figure 14 shows three normalized blade-to-blade velocity magnitude contours to represent the flow field in the second impeller at Re = 10,229 and Re = 3346, both for operations at BEP flow rates in each case. The gray zones correspond to backflow regions, indicating areas where separation occurs. The separation zone starts near the corner between the blade suction side and the impeller shroud, due to the curvature of the impeller as explained in the Figure 13. The separation zone is then spread downstream until the impeller exit. The observed pattern is similar in both cases despite the difference in Reynolds numbers, even though the backflow zones increase slightly from Re = 10,229 to Re = 3346. Figure 15 shows velocity vector plots between the diffuser exit and the impeller entrance at the midspan surface for operation at three rotational Reynolds numbers, all at the BEP flow rate for each case. Velocity in this case is calculated assuming a stationary frame of reference (V stn ) for both the diffuser and the impeller, to better illustrate the interaction between both components. As explained earlier, the abrupt change to a 90°angle at the exit of the diffuser facilitates the formation of the recirculation close to the vane trailing edge. Even though this may be designed to deliver to the next impeller a flow as free of swirl as possible, one can observe that the recirculation formed creates a low velocity zone that advances into the impeller. It can be noticed that this zone is increased as the Reynolds number decreases for the cases being analysed, since the flow rate operation in each case is reduced in accordance. This is one of the reasons why the first impeller flow field and, in accordance, its pressure rise are quite different from those at the next stages. Figure 16 shows instantaneous and time-averaged turbulence kinetic energy plots at the midspan surface in the three stages of the ESP for operations at Re = 716, 046 and Re = 10, 229. Values of turbulence kinetic energy are normalized by ( · R 2 ) 2 . When analysing the instantaneous plots in Figure 16, one can notice that turbulence generated in the impeller is propagated downstream to the diffuser. In the time-averaged plot, this effect looks smoothed inside the diffuser, throughout the whole tangential extension. Also, when Re = 716, 046 one can observe that the normalized turbulence kinetic energy is globally lower inside the first impeller, which receives a well-oriented flow from the intake pipe. When the flow advances to the first diffuser the turbulence kinetic energy value generally increases. This may be directly related to the trends observed in Figure 15, where separation zones generated in the diffuser are spread to the next stage. The overall turbulence level for Re = 10, 229 is obviously lower than for Re = 716, 046. However, the Re = 10, 229 case shows a different pattern of turbulence kinetic energy distribution inside the pump, with high levels being located mainly at the interface between the impeller and diffuser. At this region, the blade/vane interactions give rise to greater turbulence intensity levels than in other parts of the pump, despite of the large blade/vane gaps found in this type of pumps when compared, for example, with radial pumps with vanned diffusers. One also observes that, whilst the overall turbulence levels generally increase from the first stage downstream to the third stage in the Re = 716, 046 case, the opposite effect occurs for Re = 10, 229, where the viscous effect seems to gradually smooth the turbulence generated in the first stage.

Overall comments on the results
The results presented in this paper show a reasonable amount of information to help understanding the performance degradation of centrifugal pumps by the influence of the fluid viscosity. By studying a multistage ESP, it was possible to understand how flow features associated with viscous fluids affect not just a given stage, but can also influence adjacent stages. In future works, the authors intend to use proper dimensional analysis to investigate how performance degradation can be correlated between different pumps. Some aspects covered in this study are nonetheless far from being comprehensive. For example, flow regime transition seems to play a fair role in the case of oil pumping, but literature review suggests that this topic requires more research. A deeper investigation of flow field pattern for viscous pumping could help companies to design ESPs that can handle highly viscous fluids more efficiently, since it seems that many of them still rely on the inviscid flow theory to develop pump designs.

Conclusions
Numerical studies were conducted to investigate the performance degradation of an electric submersible pump operating with fluids of various viscosities. The main conclusions of this study are presented below: (1) Comparison of numerical pressure head curves with data from literature showed good general agreement, even for operation with fluids of very high viscosity. (2) Work at low rotational Reynolds number is subject to head and flow rate degradation, which reduces the range of mixed-flow pump operation. The main mechanism behind the performance degradation is the substantial increase in friction losses. (3) The flow field pattern found for the best efficiency point changes significantly with the rotational Reynolds number. An important feature to note is the appearance of separation zones inside the impeller for operation with medium-to-low Re , a fact associated with the reduction of the flow rate in the best efficiency point due to performance degradation. For operation with a very viscous fluid, however, the separation zones are generally smoothed out giving place to a well-oriented flow. (4) Analysis of the interaction between the pump parts showed that recirculation from the diffuser can influence downstream into the following impeller, as well as turbulence generated inside the impeller can be spread downstream to the diffuser. This contributes to a significant difference in performance and flow pattern from the first stage in comparison with following stages. In addition, analysis of turbulence kinetic energy contours indicated that operations at quite different Reynolds numbers are not only subject to different overall turbulence levels, but also to different turbulence distribution patterns inside the pump stages.
By using computational fluid dynamics, one gets great flexibility to simulate a number of conditions and extract useful information that may be very hard to reach by experimental methods, especially using mixed flow type ESPs. In this sense, the numerical model used proved to handle the operating conditions tested quite well, given the purposes of the present work. However, a deeper study of performance degradation under operation with highly viscous fluids should include other pump types, as this problem seems to be very sensible to the pump geometry. Another common problem in the oil and gas industry where ESPs are widely used is the presence of a gas phase in the flow, which also causes performance degradation. Among the topics of ongoing research of the present authors, to correlate performance degradation between different pumps and viscous two-phase flow pumping may be subject for future works.