Cuttings transport: On the effect of drill pipe rotation and lateral motion on the cuttings bed

Drill pipe rotation is considered a relevant factor for cuttings transport and hole cleaning. However, in the term “rotation” is often used as a moniker for the combination of plain drill pipe rotation around its own axis and more complex lateral motion, as many laboratory setups feature an unconstrained drill string. Lateral motion is generally considered to benefit the transports of cuttings due to increased bed agitation. By means of Computational Fluid Dynamics, we have investigated the effect of synchronous and asynchronous whirling drill string motion on the cuttings bed and cuttings transport for water and a more viscous, shear-thinning fluid using the Two Fluid Model in conjunction with the Kinetic Theory Of Granular Flows and closures from soil mechanics to rheologically describe granular matter. The dynamic mesh capability of ANSYS Fluent R17.2 is exploited to account for the orbital motion of the drill string. In addition, three base cases (negative eccentric, concentric, and positive eccentric) are investigated for comparison. Whirling motion helps tremendously to disperse the solids into the main flow region and hence improves the quality of cuttings transport and hole cleaning, with synchronous whirl by far outperforming asynchronous whirl due to the cumulative tangential and radial velocities. The effect is much more prominent for water than for the more viscous, shear-thinning fluid because the latter already shows a comparatively good cuttings transport performance. Moreover, in case of the more viscous, shear-thinning fluid, the positive eccentric annulus provides an even better cuttings transport capability, if comparison is made on equivalent pressure gradients. Because of the higher viscosity level, the whirling motion reduces the axial throughput, which despite the increased bed agitation results in worse performance compared to the positive eccentric case.


Introduction
In petroleum drilling, solid particles (cuttings) are generated by the drill bit which is being pushed downhole with a certain rate of penetration (ROP). The cuttings are subsequently transported by the often shear-thinning drilling fluid through the annular space (created by the drill pipe 1 in a wellbore) to the surface, as qualitatively depicted in Fig. 1.
Adequate cuttings transport is required for proper hole cleaning, i.e. the absence of a critical cuttings bed to avoid costly downtimes in drilling due to e.g. stuck pipes. The quality of solids transport depends on many factors (A. Busch et al., 2018a, Busch et al., 2019, two of which are drill pipe rotation and eccentricity. Due to the relevance of cuttings transport to the drilling industry, these have been the subject of many experimental studies (Avila et al., 2008;e.g. Han et al., 2010;Larsen, 1990;Sanchez et al., 1999;Tomren et al., 1986) over the last decades as well as numerical, or more precisely Computational Fluid Dynamics (CFD) studies (e.g. Akhshik et al., 2015;Epelle and Gerogiorgis, 2017;Heydari et al., 2017;Pang et al., 2019Pang et al., , 2018 in recent years.

Effect of parameters on cuttings transport and hole cleaning
Negative 2 eccentricity increases the accumulation of particles at the lower side of the annulus (leading to a sediment or cuttings bed) because the narrower gap results in a local reduction in fluid velocity (Bicalho et al., 2016a;Heydari et al., 2017). At the same time, pressure loss decreases because the effective cross-sectional flow area increases. This also holds for the single-phase (SP) case, regardless of fluid type or flow regime (McCann et al., 1995). Rotation on the other hand generally increases the transport of cuttings (Duan et al., 2010;Han et al., 2010), in particular in the cases off negative eccentric configurations because the tangential velocity of the rotating pipe is acting at the position of high solid volume fractions, i.e. at the cuttings bed (Bicalho et al., 2016a;Heydari et al., 2017;Xiaofeng et al., 2014). However, this effect is dependent on the particle size as small particles will be re-entrained much easier than large ones (Duan et al., 2008;Sifferman et al., 1992) as well as the annular diameter ratio, as the effect of rotation becomes much more relevant for smaller annular gaps (Peden et al., 1990). By reducing an existing cuttings bed and thereby increasing the effective flow area, drill pipe rotation leads to a decrease in pressure loss, which is different to the SP case where rotation may increase or decrease pressure losses, depending on the flow regime and fluid (Sorgun et al., 2011). For instance, pressure losses increase with rotation for turbulent flows and decrease for laminar flows of Power-Law (PL) fluids because of the shear-thinning property of the fluid (Johansen et al., 2003 et al., 1995). On the other hand, for the case of Yield-Power-Law (YPL)/Herschel-Bulkley (HB) fluids, Erge et al. (2015Erge et al. ( , 2014 observed no significant effect of rotation for turbulent flows and a pressure loss increase or decrease in the laminar and transitional regime, depending on the magnitude of inertial forces, i.e. the spatial dimensions of the annulus as well as the viscosity level of the fluid. Often, in cuttings transport studies, the drill string is assumed to rest in a fixed position, which may be either concentric or eccentric. However, this is rarely the case in wellbores (Ahmed et al., 2010;Saasen, 2014), where the drill pipe may feature complex lateral motion patterns (Gao, 2010;Leine et al., 2002;Shyu, 1989). Rotation is a necessary requirement for lateral motions. For a given rate of rotation, a variety of lateral motion patterns may be observed depending on the three-dimensional (3D) wellbore trajectory and the particular point of the wellbore one focuses on. The flexibility of the drill pipe (Xiao et al., 2003) and the buckling of the drill string (Erge et al., 2015(Erge et al., , 2014, as a consequence of the axial force applied on the string and the bit, i.e. weight on bit, determine the local eccentricity and in combination with the drill string rotation and its experienced torque (Leine et al., 2002) as well as the hydrodynamic pressure and viscous forces of the drilling fluid (Leine et al., 2002;Xiao et al., 2003) lead to a specific lateral motion.
In general, a drill pipe rotating in a wellbore with an angular velocity ωp ¼ 2π rpmp/60 may feature a variety of lateral/orbital motion patterns, which include the absence of lateral motion, i.e. pure rotation, snaking motion, where the drill pipes climbs the annular wall to a certain extent and then falls back due to gravity, irregular motion, and whirling motion, where the drill pipe rolls or slides on the surface of the outer pipe in an clockwise or anti-clockwise manner, as detailed further in section 2.1 (Gao, 2010). While specific types of lateral motion of the drill string may cause material wear and damage (Cayeux et al., 2018) as well as an increase in average pressure drop (Erge et al., 2015(Erge et al., , 2014Khatibi et al., 2018aKhatibi et al., , 2018b and pressure oscillations , it is also reasonable to expect an increased transport of cuttings because the motion of the drill pipe additionally agitates the bed and entrains particles into the bulk of the liquid flow. However, only a very limited number of research activities have focused on the specific effect of lateral drill string motion on cuttings transport.

State of the art of lateral drill string motion experiments
In experimental studies, both rotation and eccentricity have been extensively investigated and in many laboratory setups rotation and eccentricity are truly independent parameters (e.g. Duan et al., 2010;Peden et al., 1990;Sifferman et al., 1992). However, in many other laboratories (e.g. Avila et al., 2008;Khatibi et al., 2018aKhatibi et al., , 2018bSanchez et al., 1999;Sayindla et al., 2017;Ytrehus et al., 2018Ytrehus et al., , 2015 the drill string is not entirely constrained and hence lateral, orbital and/or whirling motion may occur. Thus, the eccentricity at the point of observation is time-dependent and a function of the aforementioned parameters and, unfortunately, often undisclosed. In a review conducted by Pilehvari et al. (1999), the relevance of the role of lateral drill pipe motion is first mentioned as "the manner in which the drill string behaves dynamically" and attributed to the study of Bassal (1996) in the sense that "all previous experimental studies had limitations in simulating the real dynamics of the drill pipe". At the same time, many of these results were disseminated by Sanchez et al. (1999).
In the experiments of Sanchez et al. (1999), the drill string was only constrained at its end and hence was able to move freely in the middle where flow observations were made. Sanchez et al. (1999) showed that not the pure rotation but rather the resulting orbital motion is the reason for significant improvement of hole cleaning, both during actual drilling and after drilling when cuttings generation has ceased.
Very recently, Pang et al. (2019) studied orbital drill pipe motion by means of CFD for the case of a PL fluid and showed that orbital motion increases the cuttings transport ratio (CTR) while decreasing the pressure drop (However, higher rotating speeds cause a sharp increase in pressure drop in particular when the drill pipe orbits in the opposite direction to its self-rotation). The larger the radius of the orbital motion, the better for the CTR. Orbital motion periodically stirrs up the cutting bed by causing the core zone of the axial bulk velocity following the orbital motion (though lagging behind that) and producing secondary tangential flows and eddies. Pang et al. (2019) utilized the Eulerian-Eulerian Two Fluid Model (TFM) in combination with the Kinetic Theory of Granular Flows (KTGF), though no additional closures were employed to account for frictional effects within the dense granular media, i.e. the cuttings bed. The sliding mesh method of ANSYS Fluent R14.0 was employed to realize the orbital motion of the drill pipe.
Self-rotational speed ωp and orbital speed ωw were considered to be equivalent in magnitude, i.e. ωw ¼ ωp and ωw ¼ -ωp, and investigated in the range 0-200 rpm.
Recently, Cayeux et al. (2018) investigated the special case of Fig. 1. Cuttings transport process on an annular scale.
A. Busch and S.T. Johansen synchronous whirl, where the drill string angular frequency equals the angular frequency of the whirling motion and the pipe always faces the same side towards the outer wall of the annulus, comprising the cases of SP laminar flow of Newtonian, PL, and HB fluids. Based on the methodology used by Feng et al. (2007) they inverted the annular system by considering the outer cylinder rotating around the inner and accounted for centrifugal and Coriolis forces. The whirling motion contributes to the total pressure gradient for all fluids investigated. However, in case of the HB fluid the pressure gradient of the pure rotational cases exceeds the one of the whirling cases at approximately 80 rpm for all volumetric flow rates investigated, presumably because the whirling motion avoids plug regions at higher rotational speeds (Cayeux et al., 2018). Both Vieira Neto et al. (2012) and (Bicalho et al., 2016b(Bicalho et al., , 2016a) experimentally investigated the flow of laminar xanthan gum solutions in annuli with orbital inner pipe motion. In addition, they simulated the pressure drop using the dynamic meshing (Neto et al., 2012) and sliding mesh capabilities (Bicalho et al., 2016a(Bicalho et al., , 2016b of ANSYS Fluent and obtained a good fit between experimental and numerical results. Rotation and corresponding orbital motion of the inner tube results in more uniform flow distributions in the annulus, preventing flow stagnation in the narrow-gap regions in cases of eccentric configurations. Therefore, in the case of a partially blocked annulus with eccentricity, increasing drill pipe rotation and orbital motion is considered to improve the transport of cuttings (Bicalho et al., 2016a).
For a negative eccentricity of a highly shear thinning fluid, Demiralp (2014) investigated the effect of different whirl patterns on cuttings transport. A two-way coupling between solids and fluid by means of CFD and the Discrete Element Method (DEM) and presumably (The details of the whirling motion implementation are not disclosed) the dynamic meshing capabilities of ANSYS Fluent were employed to investigate hole cleaning for different fluid superficial velocities and drill pipe rotations. Solids concentration decreases with increasing whirling speed in all flow regimes, with synchronous whirl yielding the highest solid superficial velocity.

Purpose, scope and structure of this paper
While some numerical studies investigated the impact of whirling motion on the flow hydrodynamics (Bicalho et al., 2016a;Cayeux et al., 2018;Feng et al., 2007;Neto et al., 2012), the effect of drill pipe whirl on cuttings transport has-to our knowledge-not been quantitatively investigated, with the notable exceptions of the thesis of Demiralp (2014) and the recent study of Pang et al. (2019). Experimental investigations have often featured a laterally moving drill string; however, the lateral motion is usually a consequence of the system and the controlling parameter is simply the rotational rate of the drill pipe. Unfortunately, no clear distinction is made in the literature between plain drill pipe rotation and additional lateral, orbital, or whirling motion, in particular when it comes to the interpretation and quantification of experimental results. Even recent review papers combine these under the umbrella rotation Luft, 2014a, 2014b;Ofei et al., 2015;Xiaofeng et al., 2013). In addition, we are not aware of any study where the modes of and the parameters describing orbital motion were disseminated. This hinders quantitative comparisons.
We numerically investigate the effect of two classical whirling motion cases, namely forward and backward whirl (detailed definitions are provided in the following section 2.1), on cuttings transport for water and a more viscous, shear-thinning fluid. We then compare these results to eccentric and concentric cases with plain drill string rotation.
In the following section, a description of the drill string whirl cases as well as all relevant other parameters investigated is presented, followed by a brief summary of the physical CFD model along with required closures. SP results are then presented and compared to the experiments of Khatibi et al. (2018aKhatibi et al. ( , 2018b, followed by the presentation of cuttings transport simulation results. In the subsequent discussion, we provide explanations for the observed phenomena as well as an analysis of the strength and weaknesses of our investigations. Finally, a brief summary and outlook is given.

Drill string whirl
In general, four patterns of drill pipe motion may be characterized as follows (Shyu, 1989): 1. Synchronous whirl (SW), also known as forward whirl, where the tool joint is sliding on the hole/casing wall in such a manner that it always faces the same side towards the outer wall. Consequently, the drill pipe rotation and its whirling motion show identical angular velocities, i.e. ωw ¼ ωp.
2. Asynchronous whirl (AW), also known as backward whirl, where the tool joint is rolling on the hole/casing wall without any slip. Consequently, the whirling motion occurs in opposite direction of the drill pipe rotation and the angular velocity of the whirl motion is given by ωw ¼ -dj/doωp, where dj is the diameter of the tool joints.
3. More complex whirl, where there is slip between the surfaces such that ωw ¼ csωp, where cs is different from 1 or -dj/da, not necessarily constant and may even be different for the y and z-direction such that the motion pattern becomes a Lissajou curve. 4. Any other (seemingly chaotic) motion, where e.g. the drill string does not always remain in contact with the wellbore wall (at all times) and/or where there is slip between the surfaces of tool joints and wellbore wall.
We here focus on plain whirling motion, i.e. type 1 and 2 as characterized above, because it is easiest to parametrize. A generic framework for the whirling motion is given by a 2D oscillation equation in y and z: (1) where Ei are the dimensional eccentricities, Ai the dimensional amplitudes, and ωi the angular velocities (which are here taken as ωw ¼ ωy ¼ ωz), as depicted in Fig. 2.
While SW and AW may be straightforwardly parametrized as described above, this is not so in many experimental setups because the drill string is free to move in the test section. For instance, in the experiments of Khatibi et al. (2018aKhatibi et al. ( , 2018b, the eccentricity and amplitude are functions of the drill pipe's rotation rate and the superficial fluid velocity. The drill pipe consisted of several acrylic elements connected with flexible joints, the diameter of which was slightly larger than the drill pipe (see Fig. 17 in Appendix C). One end was connected to a dual current motor while the other end was not constrained and hence free to move. Therefore, the flexible drill string arrangement was subject to lateral motion because of the enforced rotation at one end and the systems mechanical properties (compliances of individual drill string elements and flexible joints, Coulomb and viscous friction). Khatibi et al. (2018b) showed, that in their SP experiments the observed change in vertical eccentricity Ey is largely dependent on the rotation rate of the drill pipe and to a smaller extent on the Reynolds number of the flow. We curve-fitted second-order polynomials to the available dimensional data (Khatibi, 2018), the coefficients of which are presented in Table 1.
In the experiments of Khatibi et al. (2018aKhatibi et al. ( , 2018b, the horizontal amplitude Az as well as the angular velocity ω w were not determined explicitly. However, the latter may be estimated with 2π rpm=60 based on the frequency spectra analysis of pressure readings (Khatibi, 2018;Khatibi et al., 2018a), while the former may be roughly estimated by A y = 2 (Khatibi, 2018). Table 2 summarizes the SP cases investigated for rotational and lateral-motion model validation.

Test matrix
For the cuttings transport multiphase (MP) simulations, a horizontal 8.5 in wellbore section (d o ¼ 0.216 m) with a 5.0 inch drill pipe (d p ¼ 0.127 m and tool joint diameter d j ¼ 0.168 m) was assumed. Different fluids, eccentricities/whirl types, pressure gradients and drill pipe rotation rates were investigated as summarized by Table 3 in order to represent field values.
For the sake of clarity, Fig. 3 details column e y of Table 3, where the drill pipe's AW and SW motion is defined as described in section 2.1 and equations (1) and (2).
The properties of the two types of fluids investigated in this study are given in Table 4.
In all cases, following Khatibi et al. (2018aKhatibi et al. ( , 2018b, the solids were simplified as spherical particles with diameter d s ¼ 1.2 mm, density ρ s ¼ 2650 kg/m 3 and angle of internal friction α if ¼ 45 � . The solid loading was determined such that without flow and rotation, the solids bed was filling the lower clearance for the smaller eccentricity e y ¼ À 0.54, which yields α s ¼ 0.047.

Physical model
As a MP flow model, we here apply the TFM in combination with the KTGF developed by Savage (Lun et al., 1984;Savage et al., 1996;Savage and Jeffrey, 1981) handling the loose, i.e. the collisional/kinetic regime 3 (solid volume fraction α s < α s;f ¼ 0.55) and closures from soil mechanics describing the dense regime (α s > α s;f Þ of the cuttings, hereafter termed solids. The two phases are considered as interpenetrating continua and mass continuity and momentum transport equations along with closures for the fluids and solids material functions, turbulence, and the momentum exchange terms are used to model the system. The full model description is given in 5 Appendix A. As previous investigators (e.g. Epelle and Gerogiorgis, 2017;Pang et al., 2019Pang et al., , 2018Zakerian et al., 2018), we utilize the model implementation of ANSYS Fluent (ANSYS, Inc., 2016a, 2016b), a broadly used commercial CFD code.

Fig. 2.
Geometrical and kinematical system quantities, downstream orientation. Left: Geometrical dimensions of drill pipe, tool joints and outer wall of wellbore as well as drill pipe rotation around its own axis. Right: Kinematic quantities describing whirling motion of drill pipe in wellbore.

CFD setup & numerics
For the SP simulations, a variety of meshes of varying eccentricity based on Table 1 were created. Table 7 in Appendix E provides an overview of the relevant mesh parameters and Fig. 4 provides an eccentric example of the "Intermediate (MP)" case.
The dependency of the numerical solution on the mesh resolution was firstly evaluated with SP simulations (without rotation) as depicted in Fig. 19 in Appendix E. For all meshes, the r -spacing was non-uniform in order to obtain a higher resolution close to the walls. Periodic boundary conditions (BC), i.e., what leaves the domain enters the domain, were applied to either end of the annular element in order to reduce computational efforts. The length of the computational domain was chosen as L ¼ d o for the SP, and L ¼ 3d o for the MP simulations (such that any periodicity in the solution is not influenced by the BC). 4 The "Coarse (High Re)" mesh and the "Superfine (Low Re)" mesh results differ by only 1%, while the "Intermediate" mesh result differs from the "Superfine (Low Re)" mesh result by 3.8%. However, the "Intermediate" mesh does feature a much smaller first layer more suited for larger flows/pressure gradients (and thus steeper wall gradients) and it represents the experimental data best. Fig. 20 in Appendix E shows transient results for a MP case for the "Intermediate (MP)" and "Fine" meshes. In the near-steady-state time period, the difference between the two meshes is <1% for U f and about 5% for U s . The selected mesh quality "Intermediate" provided a fair compromise between sufficiently accurate results in the form of time-averaged integral quantities (Goldschmidt et al., 2004) such as superficial velocities/mass flow rates and associated computational effort and is comparable to similarly sized grids in other studies (e.g. Epelle and Gerogiorgis, 2017;Rooki et al., 2013bRooki et al., , 2013aZakerian et al., 2018).
The plain rotation of the drill pipe around its own axis was defined as a slip velocity of the inner wall, i.e., the fluid velocity at the wall is not zero but has a magnitude and direction equivalent to the wall rotational speed of the drill pipe. In case of SP simulations, we specified a mass flow rate, while in case of MP simulations a mixture pressure gradient Δp=Δx along with the solid volume fraction α s was specified.Fluent's dynamic meshing capability was employed to deform the mesh and simulate the orbital motion of the inner pipe. The latter was defined by an User-Defined-Function (UDF) which is simply an implementation of the time derivative of equations (1) and (2) and provides the velocities of the center of gravity of the orbital motion of the drill pipe to the solver. The spring-based smoothing method, were the cell edges are treated as elastic springs, was used to update the mesh every time step.
All simulations but the SP mesh dependence cases depicted in Fig. 19 in Appendix E were performed in a transient manner, as exemplarily depicted in Fig. 20 in Appendix E. The time step was 10 À 3 s to 10 À 4 s and a second-order implicit scheme was employed. The (Phase-Coupled) SIMPLE scheme (Vasquez, 2000) was used for pressure-velocity coupling. The QUICK scheme (Leonard, 1979) was used for second-order spatial discretization and the Green-Gauss node-based gradient scheme was used to evaluate all gradients. The time discretization was implicit second order. The algebraic multigrid method with the Gauss-Seidel solver and conservative under-relaxation factor settings were used to solve the system of discretized equation.
At t ¼ 0 the solids were patched into the domain and then allowed to settle over time until a partly moving-partly stationary bed is built up, resulting in an quasi-steady state of U f and U s , as illustrated by Fig. 20. Simulations were then continued for at least five orbital motions for the purpose of data sampling.
Pre-studies showed that the mesh deformation works fine for a couple of orbital motions only. After six to twelve orbital cycles, the  mesh starts to deform non-uniformly and eventually highly skewed cells lead to negative volume and divergence. We therefore simply replaced the mesh after a preset number of orbital cycles (two to four) with the initial, where the solver transfers the current solution from the old to the new mesh using interpolation schemes. Reducing the time step size did not rectify the mesh deterioration.

Results
We first present results which to some extent validate the CFD model with available experimental data. Secondly, we present results for the industrially relevant 8.5 inch wellbore section flow case where we focus on the effect of whirling motion on cuttings transport and the transitional flow regime.

Validation with single-phase experimental data
For validating the CFD model we hereafter present SP results benchmarked with respective experimental data and friction factor correlations. We use the experimental data of Khatibi et al. (2018aKhatibi et al. ( , 2018b because besides containing data for transitional flow of PL fluids, this data set also contains data for whirling drill string motion. In their experiments, Khatibi et al. (2018aKhatibi et al. ( , 2018b used water and a shear-thinning 1 g/L polyanionic cellulose solution (PAC). A broad variety of friction factor correlations for turbulent concentric and eccentric annular flow is available in the literature. For instance, for the fully eccentric annular turbulent flow of different concentrations of drag-reducing guar gum solutions, explicit friction factor correlations as a function of the generalized Reynolds number, diameter ratio, and relative roughness are available (Dosunmu and Shah, 2015;Ogugbue and Shah, 2011). Other examples are the works of Kelessidis et al. (2011) and Pilehvari and Serth (2009) for the flow of bentonite suspensions. However, to our knowledge, no dedicated friction factor correlation for PAC solutions exists. For the case of PAC, we therefore utilize the correlations of Dodge and Metzner (1959) and Irvine (1988), which are corrected for eccentricity (Haciislamoglu and Cartalos, 1994;Haciislamoglu and Langlinais, 1990), if required. Fig. 5shows CFD results for the flow of water in a concentric and fully eccentric annulus without pipe rotation.
The results are benchmarked with the aforementioned friction factor correlations from the literature. In addition, experimental data ) is depicted. The reason for choosing this particular experimental data set is that it also contains data for whirling drill string motion.
While the CFD predictions adequately fit the empirical relations for both the concentric and eccentric annular configurations, the experimental data for the eccentric case falls slightly on top of the concentric flow data.
A zoom on the low superficial velocity range of Fig. 5 (black box) is depicted in Fig. 6.
In addition, CFD results for three different whirling motion cases are depicted along with the corresponding experimental results . Note that the maximum eccentricity and the predefined whirling motion parameters for these cases are functions of the drill pipe rotational speed and the superficial fluid velocity and only valid for the experimental setup of Khatibi et al. (2018aKhatibi et al. ( , 2018b. With increasing rotational speed and the corresponding change of whirling motion, the respective pressure gradient increases. Interestingly, the aforementioned difference between CFD and experimental results is not observed for the 100 and 200 rpm cases. However, the 300 rpm cases are heavily underpredicted by the CFD simulations.
Finally, Fig. 7 provides pressure gradient vs. superficial fluid velocity results for the flow of PAC, both for the concentric and the eccentric case.
While there is good agreement between the friction factor correlations, CFD results and, in the eccentric case, the experimental data in the lower (laminar) range of superficial fluid velocities, this is not so for the higher (transitional) range. Compared to the friction factor correlations, the CFD results overpredict the transitional regime by up to 40%. In case of the eccentric configuration, CFD and experimental results show good agreement up to U f;x ¼ 1.4 m/s, where the experimental results start to deviate from the CFD results and approach the friction factor correation at higher U f;x .

Validation with multi-phase experimental data
Concerning MP flows, the physical model as presented in Appendix A and its implementation in Fluent as used in this study has been validated to a good extent by several other researchers as depicted in Fig. 8 for the case of non-whirling flow cases and based on the non-dimensional Π-space of Busch et al. (2019). Except for the Cuttings Transport Ratio (CTR) and, to some extent, the Taylor number Ta, our parameter space as given by Tables 3 and 4 is encompassed in the spaces of previous studies. The lower CTR is a consequence of our comparatively high fluid superficial velocities. Most previous studies have validated their models with cases where the drill pipe is not rotating (→ Ta ¼ 0), the reason is-from our point of view-that in many experimental studies the drill pipe is actually allowed to move freely but unfortunately this is often not clearly communicated. Therefore, high-quality experimental data suited for validation purposes is scarce.

Cuttings transport studies
While there are many ways to quantify the efficiency of cuttings transport and hole cleaning (A. Busch et al., 2018a), we here apply the CTR as the ratio of the two superficial phasic velocities (Bourgoyne et al., 1991, p. 178 where the superficial velocities are defined as where u i are the respective phasic x -velocity components and A is the cross-sectional area. 5 This CTR choice is mainly motivated by the specification of the mixture pressure gradient Δp=Δx and the solid volume fraction α s as input parameters in our numerical simulations due to the periodicity of our computational domain. Hence, the latter constitutes a fixed mass of solids and hence predetermined bed height in the absence of flow. Therefore, the superficial velocities and the pressure gradient constitute the response of the system. For clarity, it is important to realize that the results presented do compare to each other in terms of dp=dx equivalence only and not in terms of equivalence of U s f (or U s s ). While the latter is often used in the literature and is beneficial because the flow rate is known a priori, from a controls engineering point of view the former is sounder: While the volumetric fluid flow rate may be the primary variable to manipulate during operations, it is the pressure gradient which is monitored and to be kept below critical values.
For the different cases defined and depicted in Table 3 and Fig. 3, respectively, and H2O as the fluid phase, Fig. 9 shows the CTR plotted vs. drill pipe rotation and dp=dx.
In the absence of whirling motion, the CTR is highest for the positive eccentricity and lowest for the negative eccentricity. The effect of plain 5 Note that as the phases are incompressible the definition of the superficial velocities provided in equation is equivalent to the volumetric flow rates, i.e.
A. Busch and S.T. Johansen drill pipe rotation is generally highest in the case of negative eccentricity. However, for higher pressure gradients, it similarly increases the CTR for the positive eccentricity case.
The presence of whirling motion significantly increases the CTR for rotational rates faster than 60 … 100 rpm. AW leads to CTR levels between concentric and positive eccentric drill pipe arrangements, while SW is outperforming all other cases for rotational rates faster than 60 … 100 rpm.
For the different cases depicted in Fig. 3 and PAC as the fluid phase, Fig. 10 shows the CTR plotted vs. drill pipe rotation and dp= dx.
As in the Newtonian case, the CTR is highest for SW, given that rotational rates larger than 60 … 100 rpm are maintained. As opposed to the Newtonian case, the AW case falls between the negative and concentric cases for the entire range of rotational rates.
The effect of plain drill pipe rotation is largest for the negative eccentric case and virtually non-existent for the positive eccentricity. For the concentric case, the CTR jumps from one level to another between 30 and 60 rpm.
In addition to the CTR ¼ f(rpm, dp=dx), we provide the results in the form ROP ¼ f(rpm, dp=dx), where ROP is related to the superficial solid velocity as the nominator of the CTR as follows: In a real field scenario at steady-state (with respect to all input parameters such as U s f and ROP), the superficial velocity of the solids U s s is determined by the ROP, the bit diameter D o , and the rock porosity α r as a consequence of mass conservation.  Khatibi et al. (2018aKhatibi et al. ( , 2018b. In addition, CFD results for a concentric annulus are depicted. Empirical correlations for both the concentric and eccentric-corrected for eccentricity (Haciislamoglu and Cartalos, 1994;Haciislamoglu and Langlinais, 1990)-case are plotted with dashed lines. The black box highlights the area depicted in Fig. 6.  . In addition, CFD results for a concentric annulus are depicted. Empirical PL correlations for both the concentric and eccentric-corrected for eccentricity (Haciislamoglu and Cartalos, 1994;Haciislamoglu and Langlinais, 1990)-case are plotted with dashed/dotted lines.  Fig. 9. Top: Absolute CTR vs. drill pipe rotation rate and pressure gradient for H2O and the systems defined in Fig. 3. Bottom: Relative change of the CTR based on the concentric system e0.

A. Busch and S.T. Johansen
While we find the same qualitative results for the case of H2O (Fig. 11), this is not so for the case of PAC (Fig. 12), where in contrast to the CTR ¼ f(rpm, dp=dx) presentation of Fig. 10 the positive eccentricity yields the highest ROP for the entire range of rotational rates considered and hence performs best in terms of hole cleaning.
To further illustrate the effect of the varied parameters on the results, we additionally depict the fluid superficial velocities for H2O and PAC in Fig. 13 and Fig. 14, respectively. As expected, due to the different viscosity magnitudes, the fluid throughputs of PAC consistently fall under the levels of H2O for a given dp=dx.
For the concentric and positive eccentric cases, the effect of rotation is a bit more pronounced at lower dp=dx, and correspondingly fluid superficial velocities. For the negative eccentric case, as opposed to the H2O, a plateau exists for rotation rates >100 rpm in case of the PAC solution.
For any given pressure gradient, orbital motion results in a significantly reduced throughput for the entire range of rotational rates considered.

Validation with single-phase experimental data
For the SP H2O base case (Re ¼ 4900, e ¼ À 0.95, 0 rpm) of Khatibi et al. (2018a), the CFD results do fairly well fit the empirical pressure drop correlation (Blasius, 1912;Haaland, 1983) with the eccentricity correction (Haciislamoglu and Cartalos, 1994) applied. However, the experimental results (Khatibi et al., 2018a) exceed the CFD results by Fig. 10. Top: Absolute CTR vs. drill pipe rotation rate and pressure gradient for PAC and the systems defined in Fig. 3. Bottom: Relative change of the CTR based on the concentric system e0.

A. Busch and S.T. Johansen
50% and do coincide with the fully concentric CFD results, which also match well with the empirical correlation (Blasius, 1912;Haaland, 1983). The significantly larger pressure drops found in the experiments may be a result of different factors. First, the computational domain assumes periodicity along x, which is not necessarily the case because of development length effects and a likely skewed whirling motion of the drill string.
The predefined whirling motion of the drill string is not necessarily describing the motion of the drill string for every x -location of the annulus. In the CFD model, we assume that the axis of the string and the outer pipe are parallel. However, due to the compliance of the drill pipe material as well as the flexible joints and the concentric fixation of the drill string at the motor end it is very likely that the drill string in the experiments features more complex whirling motion, which additionally varies in the streamwise direction, i.e. is skewed along x. Closer to the motor end it will naturally feature a more concentric and less whirling motion while further away of the motor it may move more freely and hence feature more complex elliptic motion patterns as indicated in Fig. 1 of Khatibi et al. (2018a). This is corroborated by the geometrical constraints introduced by the flexible drill string section joints with an outer diameter d j ¼ 0.031 m (Khatibi, 2018) which yields a dimensionless eccentricity e ¼ À 0.6, whereas the factual eccentricity for the no flow/no rotation situation as reported by Khatibi et al. (2018b) was e ¼ À 0.94 (see tabulated data in Table 6/Appendix C). In addition, the Fig. 11. Top: Absolute ROP vs. drill pipe rotation rate and pressure gradient for H2O and the systems defined in Fig. 3. Bottom: Relative change of the ROP based on the concentric system e0.

A. Busch and S.T. Johansen
parameters characterizing the whirling motion of the drill pipe, namely the y -and z -amplitude and the frequency, were not precisely measured but rather estimated based on the obtained experimental data. The data characterizing the vertical eccentricity and amplitude of the drill string does not cover the entire parameter space and was obtained by graphical analysis of the PIV pictures. Both the horizontal amplitude and the frequency of the whirling motion where simply estimated based on visual observations rather than directly measured.
Furthermore, hydrodynamic entrance effects may be of relevance. For the laminar flow of water in a concentric annulus with d i = d o > 0.5, the development length may be estimated with :631 1:6 þ ð0:0442ReÞ 1:6 Þ 1=1:6 (Poole, 2010), which gives 1.71 m for the Re ¼ 4900 case and 2.69 m for the Re ¼ 7700 case. In contrast, in case of turbulent flow, the development length is much shorter and may be estimated with 4:4d h Re 1=6 (Çengel and Cimbala, 2006), which yields only 0.27 m and 0.29 m, respectively. However, Lien et al. (2004) recommend 150d h =2, which yields 1.25 m.
In any case, here we are dealing with transitional flow, which is intermittent in the sense that both laminarization and development of turbulence are competing. In the experiments of Khatibi et al. (2018aKhatibi et al. ( , 2018b, the distance from the beginning of the annular section to the first pressure transducer is 1.4 m, and 1.52 m from the first to the second pressure transducer. Moreover, 0.3 m prior to the first pressure transducer, the first flexible joint with an outer diameter d j ¼ 0.031 m (Khatibi, 2018) significantly reduces the cross-sectional flow area and hence introduces a flow disturbance. Thus, the flow may still have been of developing nature in the section where pressure measurements were taken. Finally, in the case of PAC, the discrepancies between CFD and experimental results as well as the friction factor correlations may be attributed to the viscoelastic and/or drag-reducing capabilities of the PAC solutions utilized (Alexander Busch et al., 2018b). The Generalized Newtonian Fluid framework with the PL material function does neither account for normal stress differences nor for elongational viscosity, both of which affect flows in eccentric annuli with the latter also affecting rotational flows (Escudier et al., 2002). Moreover, while the employed k-ω SST model is versatile regarding yþ values it does neither consider n PL -dependent damping functions (e.g. Malin, 1997) nor non-Newtonian wall functions (e.g. Johansen and Mo, 2015) in the y þ < 1 and y þ > 30 regimes, respectively.

Cuttings transportplain drill pipe rotation
Focusing on the plain drill pipe rotation cases first, the positive eccentric case generally leads to much better hole cleaning than the eccentric case, in line with other studies (Bicalho et al., 2016a;Heydari et al., 2017). The more clearance between the drill pipe and the cuttings bed, the better the hole cleaning because of the higher fluid velocities below the drill pipe and on top of the sediment bed. Consequently, the shear stress acting on the bed is much higher for a positive eccentric drill pipe than for a negative eccentric one, hence the better CTR and ROP values. This is physically sound and in accordance with the often-stated order when it comes to the relevance of individual cuttings transport parameters: Volumetric fluid flow rate is typically considered the most important parameter, simply because it is just the axial flow components which transports solids. Drill pipe rotation is an additional contributing factor which depends on the flow regime, fluid rheological parameters, and, as shown, eccentricity. It is important to note that the results were obtained for a total solid volume fraction α s ¼ 0.047. Larger values will lead to a higher cuttings bed in the computational domain, which will  The effect of plain drill pipe rotation is highest for the negative eccentricity because the corresponding tangential velocities of the fluid and subsequently the solid phase act directly on the bed, agitate the bed and disperse solids into regions of higher fluid velocity, where they are easily transported downstream. This effect is less prominent in the other configurations. However, it is also less relevant for hole cleaning in the other configurations because of the aforementioned role of the locally higher fluid streamwise velocities and the higher shear stresses acting on the bed.
For a shear-thinning fluid in laminar flow, rotation may generally reduce the pressure gradient, i.e. increase the throughput for a given dp= dx. In contrast, in turbulent flow, rotation may increase the pressure gradient, i.e. decrease the throughput for a given dp= dx, as it increases the degree of turbulence. For the investigated shear-thinning PAC solution, this effect is very small, presumably because shear-thinning and turbulence generation are counteracting each other as the flow is in fact transitional rather than fully turbulent (see Fig. 18 in Appendix D). In addition, the presence of solids certainly overshadows this SP effect.

Cuttings transportwhirling motion
Assuming a concentric drill pipe as the base case, the presence of whirling motion generally increases both the CTR and the ROP. The same mechanism as for plain drill pipe rotation applies: The additional whirling motion leads to an increase in tangential and here additionally radial (with respect to the streamwise flow direction) velocities which help to agitate the sediment bed and disperse cuttings into the main flow regions and thus enhance cuttings transport. Note that in case of the shear-thinning PAC, the SW and AS ROP is less than for the concentric case and only for higher Δp=Δx SW outperforms the concentric arrangement. The orbital drill pipe motion leads to two counteracting effects: (1) A reduction of the viscosity due to the additional applied shear and (2) the agitation of the bed due to the increase of turbulence and increased tangential/radial velocities, the.second one becoming dominant for higher Δp=Δx.
The SW significantly outperforms the AW because in case of SW both the tangential velocity induced by the drill pipe motion around its own axis and the radial velocity induced by the drill pipe orbital motion act in the same direction and are therefore additive. In case of AW, the plain drill pipe rotation is in the opposite direction of the drill pipe orbital motion and the respective velocities to some extent counteract each other.
In case of the PAC, the solids transport is highest for the positive eccentric case. For this geometrical arrangement, the rotation of the drill pipe around its own axis leads only to a minor shift of the cuttings bed towards the side of the annulus (Fig. 16 (top) in 5 Appendix A), with the majority of the helical flow pattern occurring on top of the main flow field. This is not so for the SW case (Fig. 16 (bottom) in 5 Appendix A), where the orbital motion of the drill pipe leads to a circumferentially alternating helical flow pattern affecting the entire volume of the annulus. While this is generally considered a positive feature in the sense of bed agitation, it also leads to a reduction in fluid throughput because much more fluid obeys the induced helical motion (Fig. 14). This is also the case for H2O SW (Fig. 13), where the difference lies in the magnitude of velocities and viscosities associated with the two fluids. In case of H2O SW, the circumferentially alternating helical flow pattern leads to much less tangential flow (Fig. 15 in 5 Appendix A) and hence allows for more throughput. This indicates that (synchronous) whirling motion is most effective in enhancing cuttings transport in presence of low viscosity fluids and larger rotation rates. However, more comprehensive data is required in order to adequately assess the coupled effect of whirling motion and fluid rheological properties.
An operational challenge currently discussed in the drilling industry is drilling with an ROP of 60 m/h and a drill pipe rotation rate of 60 rpm (Iversen and Islam, 2018). Figs. 11 and 12 suggest that this may be achieved by ensuring a positive eccentricity or SW state of the drill string in the wellbore in combination with a dp=dx > 500 Pa/m for H2O and dp=dx > 250 Pa/m for PAC. However, the quantitative results presented in Figs. 11 and 12 may not simply be applied to field scenarios because the modelling framework used is, for the multiphase part, neither validated nor tuned with experimental data.

Strength and weaknesses
The presented model and computational approach is a comparatively simple tool to analyze the effect of orbital drill pipe motion. As applied in this study, it allows for quantification of the effect of whirling motion on cuttings transport and qualitatively confirms the conclusion of Sanchez et al. (1999) that (if compared to an negative eccentric drill string arrangement) the orbital motion of the drill pipe is the major reason for significant improvement of cuttings transport.
While the model and code implementation as utilized in this study has been validated to a good extent by several other researchers (e.g. Amanna and Khorsand Movaghar, 2016;Epelle and Gerogiorgis, 2017;Kamyab and Rasouli, 2016;Pang et al., 2019Pang et al., , 2018 for the case of non-whirling flow cases (see Fig. 8), further validation work is required for the whirling cases. However, this requires experimental data where the kinematics of the orbital drill pipe motions are clearly quantified. i. e., a precise description of the drill string orbital motion is provided. If this is not so, any then required estimate of motion-relevant parameters likely leads to bad model predictions as the comparison of SP simulations and experimental results of Khatibi et al. (2018aKhatibi et al. ( , 2018b) has shown.
The design space must be analyzed more comprehensively. In terms of fluid rheological properties, more viscous fluids have to be investigated as well as the role of a potential yield stress in the presence of whirling motion. In addition, solid volume fractions, solid particle diameters, the pipe-hole-diameter combination and inclination, which is a critical parameter as it defines the potential for avalanches, is to be varied in order to obtain a broader quantitative picture of the relevance of whirling motion. Given enough data, one may then also transform the data easily into the more common (α s , dp=dx) ¼ f(U f , …) framework. More complex motion patterns need to be studied as we have just focused on easy-to-parametrize forward and backward whirl. Any slip between the drill pipe collars and the wellbore wall as well as detachment may occur in the wellbore. Furthermore, the presence of the cuttings bed will likely change the circular or elliptical orbital motion patterns typically utilized in the industry (due to their simple mathematical description). However, a bidirectional coupling of moving drill pipe structure and flow of fluid and solid phases is not reasonably possible on the annular scale as the structural deformation and associated non-flow forces depend on information up-and downstream of the annular domain under investigation.
While the applied GNF framework is the state-of-the-art in cuttings transport research, it does not account for potentially relevant physics such as thixotropy and viscoelasticity. Laboratory fluids such as CMC and PAC are known to act thixotropic and viscoelastic (Alexander Busch et al., 2018b), and the viscoelastic properties of drilling fluid systems as used in the field may lead to sediment bed cohesion (Werner, 2018) that is not captured by the GNF framework.
For higher superficial velocities, the model does not replicate pressure drop quantitatively well. Two effects come into play: (1) Too high values of the fluids viscosity are to be expected due to the utilized turbulence modelling approach, which will reduce particle settling and the mass flow rate for a given pressure gradient. As briefly mentioned in 5A.4, the RANS framework of commercial solvers in general and Fluent specifically does not account for the viscosity as a varying quantity. Generally, in RANS turbulence models, the rate of strain is defined as the symmetric part of the mean velocity field gradient. This neglects any additional variation due to the velocity fluctuations, which will lead to an underestimation of the rate of strain magnitude and thus A. Busch and S.T. Johansen overestimation of the fluid viscosity. Note that the same holds for other two-equation models such as the often-utilized k -ε -model. (2) The high-Re approach taken is based on Newtonian wall functions which will likely produce incorrect estimates of the respective near-wall quantities in the non-Newtonian case.
The transport rate of solids through the domain may be overestimated by an unknown extent. The Eulerian-Eulerian method employed in this study, i.e. the Two Fluid Model (TFM) continuum approach with the Kinetic Theory of Granular Flow (KTGF) and additional closures to handle the dense granular regime, as implemented in Fluent R17.2, does not produce angle of repose satisfying conditions in the absence of flow under all conditions (Busch and Johansen, 2018). Even for a horizontal bed under the sole influence of gravity, the top-layer always remains in a liquid-like state regardless of flow time. The KTFG dynamics act in a checkerboard-like manner ensuring a very low viscous solid phase, which in the concept of the TFM continuum approach should feature high solid viscosity levels representing the non-flowing sediment bed. According to these observations, it is expected that the solids bed will behave as "fluidized" in the simulations and that the solids flux may be overestimated. An alternative and with respect to the above mentioned overestimation of the solids transport rate better-suited modelling approach is the CFD-DEM framework (e.g. Akhshik et al., 2015;Zhang et al., 2016), which if combined with the periodic BC approach as used in this study may also allow reasonable computation times.
Another issue with the TFM-KTGF-SM approach is the potential violation of its inherent continuum assumption for specific combinations of system sizes and particle sizes (Goldschmidt et al., 2004). For instance, the smallest mesh size as a consequence of the dynamic meshing technique employed is approx. 1.04 mm, which is in the order of the particle size d p . This may lead to an error in regions of high fluid velocity or pressure gradients since interaction forces between phases will be simply computed based on the respective cell values. However, the same applies for Lagrangian methods. Alternatives, such as Fluent's Macroscopic Particle Model (Agrawal et al., 2004), are not fit for purpose due to the associated computational effort for the systems under consideration here.
The employed mesh moving feature of Fluent R17.2 led to severe mesh deformation with time. Depending on the orbital frequency, the mesh had to be replaced (i.e. the current solution is interpolated to a new mesh) after several seconds of flow time in order to avoid grid deterioration. The mesh motion feature of Fluent is not meant to be used for high rotational mesh deformation, therefore a sliding mesh approach as employed by Bicalho et al. (2016aBicalho et al. ( , 2016b) may be the better choice.
The accuracy of the results may be increased by refining the mesh. This may be achieved at no additional computational costs by significantly shortening the domain since no streamwise development of any quantity was observed.

Conclusions & outlook
We have numerically investigated the role of whirling drill string motion on cuttings transport by means of CFD and a dynamic mesh technique. The essential findings are: � In case of a negative eccentric annulus, whirling motion helps tremendously to disperse the solids into the main flow region and hence improves the quality of cuttings transport and hole cleaning. The effect is much more relevant for water than for the investigated more viscous, shear-thinning fluid because the latter already shows a good cuttings transport performance. � Synchronous whirl is much better suited to agitate the bed and disperse cuttings than asynchronous whirl because the tangential and radial velocities add to each other. � For the investigated parameter values, the positive eccentric annulus provides an even better cuttings transport capability for PAC being the carrier fluid. Whirling motion reduces the axial throughput, which despite the increased bed agitation results in worse performance compared to the positive eccentric case. � The classical view of rotation being a relevant parameter for cuttings transports needs to be detailed: The cuttings transport research community needs to distinguish between plain drill pipe rotation around its own axis and rotation involving different types of whirling motion or more complex lateral motion patterns. Experimentalists are advised to carefully design their laboratory setups such that occurring whirling motion can be quantified. � More research is required to explore the entire industry-relevant design space, i.e. other numerical values of the solid volume fraction, other fluids, inclination, and orbital motions. In addition, the laminar flow regime needs to be addressed. However, a bi-partisan approach is need where experimental work is conducted in order to validate and improve the simulation work.

Appendix A1. The Cauchy equations of motion for a two-phase flow
In the TFM framework, the fluid (index f) as well as the solid (index s), phase are described as interpenetrating continua. Both fluid and solid are considered isothermal and incompressible. 6 For an arbitrary volume element V i , the phase volume fractions α i must therefore sum to one.
Mass conservation is given by where the index i 2 ff; sg and ρ i and u i denote the intrinsic volume averages of density and velocity, respectively.
Both phases obey a general form of the Cauchy momentum transport equation, which for the fluid and solid phase respectively reads where τ i is the phasic deviatoric stress tensor comprising some constitutive equation, here a compressible Generalized Newtonian Fluid (GNF) and phase-dependent material functions for the shear and bulk viscosities, η i and κ i , where D i is the symmetric part of the fluid or solid velocity gradient (also known as the rate of deformation tensor, or alternatively the rate of strain tensor) and the shear rate _ γ i is the magnitude of the rate of deformation tensor D i , _ γ i ¼ ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi The closures for the granular phase are provided in section A.3 and the rheological closures of the fluid are provided in section A.2. As we are not solving these balance equations to the smallest length scales of the flow, the phenomenon of turbulence is to be modeled. Performing Reynolds averaging (Reynolds, 1895) of the instantaneous balance equations for mass and momentum, a so-called Reynolds stress term r⋅ ðÀ ρu ' u ' Þ arises in the now ensemble-averaged momentum conservation equation. The Reynolds or turbulent stress tensor τ i;t ¼ À ρ i u ' i u ' i is usually assumed symmetric and may be modeled by applying the Boussinesque (1877) hypothesis, also known as the gradient diffusion hypothesis, to relate the Reynolds stresses to the mean velocity gradients and the turbulent viscosity in the form of The employed closures for the turbulent (or eddy) viscosity μ i;t used in the constitutive equation for the turbulent stress tensor τ i;t are further detailed in section A.4.
The last terms in equations (8) and (9) are representing the momentum transfer of one phase to the other, where the force sum is to be taken over all particles in the volume V. We here only consider the drag force f D , which is typically modeled based on the relative velocity as To model the interphase exchange coefficient K, we apply the formulation of Gidaspow et al. (1992), which is a combination of the Wen and Yu (1966) model and the Ergun (1952)  We are here concerned with either Newtonian (constant viscosity, e.g. η f ¼ 0.001002 Pa s for water) or purely shear-thinning fluids which obey the GNF constitutive equation (10) with κ f ¼ 0 assuming incompressibility. In case of shear-thinning fluids the most simple formulation of the viscosity η f accounting for shear-thinning behavior is the Ostwald (1925) material function, also known as power law (PL), Drilling fluids may feature a yield stress and are therefore typically described with the Herschel and Bulkley (1926) material function, also known as Yield Power Law (YPL). However, we here limit ourselves to PL fluids as the experimental data used for SP validation is based on a PL fluid.

Appendix A.3. Rheological closures of the solid phase
If the TFM-KTGF framework is used to describe dense granular flows, the entire solid stress tensor, namely equation (10) with index s and including the solid pressure p s , is given by the sum of collisional/kinetic and frictional components (Savage, 1983) Even though the general stencil is that of a compressible Newtonian fluid, namely equation (10), the rheological properties of the solid phase given by the respective material functions as summarized in Table 5 are highly non-linear as they depend on a variety of variables.

Table 5
Overview of solid phase state equations and material functions used to model the kinetic/collisional (index k=c) and frictional (index f ) regimes.
Schaeffer (1987) Bulk viscosity n/a n/a For instance, for vanishing shear rates in the frictional regime, a A. Busch and S.T. Johansen Bingham-type flow behavior is obtained due to the yield feature inherent in equation (26). The kinetic/collisional closures given in Table 5 are functions of the granular temperature Θ s as a measure for the degree of random particle motion, for which the transport equation reads (Ding and Gidaspow, 1990) (27) where k Θs is the granular conductivity (e.g. Syamlal et al., 1993) and the granular temperature Θ s , a measure for the granular fluctuations due to individual particle collisions, is defined as Here, u s;i ' is the i -th fluctuating component of the solids velocity in the Cartesian coordinate system and the bracket represents an ensemble average of the fluctuating velocities of all particles within a finite volume and time period (Ding and Gidaspow, 1990).
The partial differential equation (21) can be simplified to an algebraic equation by neglecting the convection and diffusion termsan often used assumption in dense, slow moving fluidized beds, where the local generation and dissipation of granular temperature far outweigh the transport by convection and diffusion. The two final terms in equation (21) are the collisional dissipation of energy (Lun et al., 1984) and the interphase exchange between the particle fluctuations and the liquid phase (Gidaspow et al., 1992).
In equations (21)-(24), e ss ¼ 0.9 is the coefficient of restitution for particle collisions and is the radial distribution function accounting for the probability of particle collisions, which has been used frequently in the history of granular flows (Bagnold, 1954;Lun et al., 1984;Ogawa et al., 1980;Sinclair and Jackson, 1989) in the form presented in equation (23).

Appendix A.4. Turbulence closures
Concerning the fluid phase, we here use the Shear Stress Transport (SST) k -ω model (Menter, 1994(Menter, , 1993, because of its suitability for swirling flows, the possibility to either integrate it to the laminar sublayer or apply wall functions, and because it correctly collapses to the laminar solution in case of laminar flows.
Dropping the fluid index f, the turbulent viscosity is defined as where l is a limiter coefficient ensuring that overprediction of the turbulent viscosity is avoided and therefore enabling the SST k -ω model to better predict the onset and amount of flow separation from smooth surfaces.
The two transport equations for the turbulent kinetic energy k and the specific dissipation ω are ∂ ∂t and ∂ ∂t where σ i are the respective turbulent Prandtl numbers, G i are respective production terms, Y i are respective dissipation terms and D ω is a cross diffusion term, which arises in equation (26) as a consequence of the blending of the standard k -ω model and the standard k -ε model (converted to a k ω formulation). For further details as well as all relevant closures of the model, the reader is referred to Menter (1994Menter ( , 1993. The solid phase turbulence is also described with a turbulent viscosity, i.e. equation (13). However, the turbulence quantities of the solid phase are A. Busch and S.T. Johansen obtained from the fluid phase by applying the Tchen theory of dispersion of discrete particles by homogenous turbulence.as given by Simonin and Viollet (1990).
To our knowledge, no non-Newtonian modifications-for instance in the form of damping functions as for the k -ε model (Malin, 1997) or in the form of additional closures for averaged apparent viscosity and turbulent cross-correlations with fluctuating viscosity as a consequence of Reynolds-averaging of the non-constant viscosity (Gavrilov and Rudyak, 2016)-have so far been developed for the k -ω family of models. Hence, we employ the SST k -ω model as implemented in Fluent, and some inaccuracy is expected in the case of non-Newtonian liquids. Fig. 15. α s and u s fields for eþ (top) and SW (bottom, lower drill pipe position depicted, full video available at https://youtu.be/vw4LUL3dF-c), H2O, dp= dx ¼ À 500 Pa/m and 130 rpm.
A. Busch and S.T. Johansen Appendix C. Experimental data Fig. 17. Annular test section of .

Table 6
Eccentricity E y and amplitude A y as observed by  and courtesy of Khatibi (2018). 0.5354 300.0000 n/a n/a n/a n/a n/a 9.400E-04 1.2275 0.0000 À 0.0061 0.0000 À 0.0061 À 0.82 À 0.82 9.400E-04 1.2275 100.0000 À 0.0056 0.0035 À 0.0020 À 0.74 À 0.27 9.400E-04 1.2275 200.0000 n/a n/a n/a n/a n/a 9.400E-04 1.2275 300.0000 n/a n/a n/a n/a n/a Fig. 18. Metzner and Reed (1955) Reynolds number vs. drill pipe rotation rate and pressure gradient for H2O (top) and PAC (bottom). See Table 3 and Fig. 3 for the test matrix and the system definition, respectively.   Fig. 19. Mesh dependency of SP water flow (no rotation) for the meshes defined in Table 7 and a pressure difference of Δp=Δx ¼ 30 Pa/m. When comparing the CFD results to the Blasius friction factor correlation, the difference is in the order of À 3 … -6% (dashed brown curve). When comparing the simulation results to the low Reynolds number superfine mesh where the wall layer is fully resolved, the difference is in the order of À 1 … 4% (solid brown curve). (For interpretation of the references to colour in this figure legend, the reader is referred to the Web version of this article.).