The Effect of Thixotropy on Pressure Losses in a Pipe

Abstract

Drilling fluids are designed to be shear-thinning for limiting pressure losses when subjected to high bulk velocities and yet be sufficiently viscous to transport solid material under low bulk velocity conditions. They also form a gel when left at rest, to keep weighting materials and drill-cuttings in suspension. Because of this design, they also have a thixotropic behavior. As the shear history influences the shear properties of thixotropic fluids, the pressure losses experienced in a tube, after a change in diameter, are influenced over a much longer distance than just what would be expected from solely entrance effects. In this paper, we consider several rheological behaviors that are relevant for characterizing drilling fluids: Collins–Graves, Herschel–Bulkley, Robertson–Stiff, Heinz–Casson, Carreau and Quemada. We develop a generic solution for modelling the viscous pressure gradient in a circular pipe under the influence of thixotropic effects and we apply this model to configurations with change in diameters. It is found that the choice of a rheological behavior should be guided by the actual response of the fluid, especially in a turbulent flow regime, and not chosen a priori. Furthermore, thixotropy may influence pressure gradients over long distances when there are changes of diameter in a hydraulic circuit. This fact is important to consider when designing pipe rheometers.

1. Introduction

Drilling fluids play multiple and important roles in drilling operations, such as maintaining an adequate pressure in the open hole section of the borehole to avoid formation fluid influx or wellbore instabilities, transporting drill-cuttings to the surface, and cooling down the downhole equipment. Drilling fluids need to be sufficiently viscous to transport drill-cuttings and yet the pressure drops caused by viscous flow should be limited to avoid excessive pump pressures. As drilling fluids are pumped into a drill-string and returned through an annulus offering a larger cross-sectional area than the drill-pipes, they are designed to be shear-thinning, therefore limiting viscous pressure drop at high bulk velocity while presenting a high effective viscosity at low bulk velocities. The downhole hydrostatic pressure is related to the density of the drilling fluid, which often results from the addition of weighting solid particles in the fluid mix. Drilling fluids are therefore designed to gel when left at rest in order to keep both weighting material and drill-cuttings in suspension, thus exhibiting a viscoelastic behavior. As a result of these design principles, their rheological properties depend on the shear history, which means that in practice they are thixotropic.

The estimation of viscous pressure losses associated with the circulation of drilling fluids has been calculated for many years, but mostly without considering the impact of thixotropy [1]. Yet under pump acceleration, thixotropy must play a role on the time evolution of pressure along the hydraulic circuit [2,3]. Additionally, any change in diameter along the hydraulic circuit influences the shear history of the fluid and therefore the resulting viscous pressure losses [4,5].

2. Viscous Flow in a Pipe

The flow of a fluid in a cylindrical pipe is axisymmetric. We will consider the force balance, in the axial direction and in laminar flow conditions, on a small shell element of radius $r$ (dimension [L] (m)), width $dr$ (dimension [L] (m)) and length $ds$ (dimension [L] (m)) (see Figure 1):

$$\frac{d\left(\tau r\right)}{dr}=r\frac{dp}{ds},$$

(1)

where $\tau$ is the shear stress (dimension [ML⁻¹T⁻²] (Pa)) and $p$ is the pressure (dimension [ML⁻¹T⁻²] (Pa)).

Since the viscous friction component of pressure does not depend on the radial position in a cross-section, we can integrate Equation (1):

$$\tau \left(r\right)=\frac{r}{2}\frac{dp}{ds}+\frac{C}{r},$$

(2)

where $C$ is an integration constant. At $r=0$, $\frac{C}{r}$ is a finite value only when $C=0$, and therefore the shear stress at a radial position is simply:

$$\tau \left(r\right)=\frac{r}{2}\frac{dp}{ds},$$

(3)

and more particularly, at the pipe wall, we have:

$${\tau }_w=\frac{d_w}{4}\frac{dp}{ds},$$

(4)

where ${\tau }_w$ is the shear stress at the wall and $d_w$ is the internal pipe diameter.

Equation (3) means that the shear stress is a linear function of the radial position. However, certain fluids may exhibit a yield stress (here noted ${\tau }_{\gamma }$), i.e., a minimum stress for which viscous flow cannot happen. For such fluids, the shear stress remains constant and equal to the yield stress when the radial position is smaller than a certain value, here noted $r_p$. As it is possible to obtain different values of the yield stress as a function of the measurement method, it is debated whether the yield stress really exists [6]. For the sake of this paper, we will not take sides, and consider both cases, i.e., with and without yield stress. The viscous shear stress is in turn a function of the shear rate, $\dot{\gamma }$ (dimension [T⁻¹] (s⁻¹)), that characterizes the slippage between two adjacent shells. The shear rate at the wall is noted ${\dot{\gamma }}_w$, and the shear rate at the center of the pipe or in the plugged region, i.e., $r\le r_p,$ is necessarily zero. The local fluid velocity, $v$ (dimension [LT⁻¹] (m/s)), is related to the shear rate by: $\dot{\gamma }=\frac{dv}{dr}$. The fluid velocity at the wall is equal to zero as a result of the no slip at the wall hypothesis, while the fluid velocity is maximum at the centerline or in the non-shear region. Figure 2 shows a schematic view of the shear stress, shear rate and fluid velocity in a cross-section of the tube.

Finally, the volumetric flowrate ($Q$, dimension [L³T⁻¹] (m³/s)) is related to the fluid velocity field by:

$$Q={{\displaystyle \int }}_0^{\raisebox{1ex}{$d_w$}\!\left/ \!\raisebox{-1ex}{$2$}\right.}2\pi v\left(r\right)rdr.$$

(5)

In laminar flow, as everything else is imposed by the pipe dimension and physical relationships (e.g., relation between flowrate and velocity field, relation between shear rate and velocity, shear stress as a function of radial distance), the volumetric flowrate and the pressure gradient only depend on the relationship between the shear stress and the shear rate.

3. Pressure Gradients of Non-Thixotropic Fluid in a Pipe

The simplest relationship between shear stress and shear rate is the Newtonian rheological behavior:

$$\tau ={\mu }_N\dot{\gamma },$$

(6)

where ${\mu }_N$ is the viscosity (dimension [ML⁻¹T⁻¹] (Pa.s)). All other relationships that are not a direct proportionality between the shear stress and the shear rate are called non-Newtonian.

The relation between the shear stress and the shear rate can be measured with a rheometer. For instance, in a concentric rheometer configuration, the rotational speed is converted to a shear rate at the wall and the shear stress at the wall is deduced from the measured torque. However, the relationship between rotational speed and shear rate at the wall is constant only for Newtonian fluids. For non-Newtonian fluids, the relationship is more complex and requires the utilization of all the measurements of shear rate and shear stresses in order to be estimated [7]. Similarly, the ratio of the torque to the shear stress is constant only for a Newtonian fluid. For non-Newtonian ones, corrections for end-effects shall be applied [8].

In practice, drilling fluids are non-Newtonian, and several rheological behaviors are candidates to describe their viscous properties.

3.1. Non-Newtonian Rheological Behaviors

A first simple extension of the Newtonian rheological behavior is the Bingham plastic model:

$$\tau ={\tau }_{\gamma }+{\mu }_p\dot{\gamma },$$

(7)

where ${\tau }_{\gamma }$ is the yield point [ML⁻¹T⁻²] (Pa) and ${\mu }_p$ is the plastic viscosity with the dimension [ML⁻¹T⁻¹] (Pa.s). The effective viscosity can be defined as: ${\mu }_{eff}=\frac{\tau }{\dot{\gamma }}$ (dimension [ML⁻¹T⁻¹] (Pa.s)). For the Bingham plastic model, it is noticeable that the effective viscosity tends to infinity when the shear rate tends to zero.

To eliminate the problem associated with an effective viscosity tending to infinity at zero shear rate, the Collins–Graves model utilizes a multiplicative term that ensures a finite effective viscosity at zero shear rate [9]:

$$\tau =\left({\tau }_{\gamma }+{\mu }_p\dot{\gamma }\right)(1-e^{-t_{CG}\dot{\gamma }}),$$

(8)

where $t_{CG}$ is a time constant [T] (s).

Another simple generalization of the Newtonian rheological behavior is the power law model:

$$\tau =K{\dot{\gamma }}^n,$$

(9)

where $K$ is the consistency index (dimension [ML⁻¹Tⁿ⁻²] (Pa.sⁿ) and $n$ is the flow behavior index (dimensionless). It is somewhat disturbing that the physical quantity of the consistency index depends on $n$, as it may lead to believe that consistency indices cannot be compared with each other. This concern can easily be removed by normalizing the shear rate. Indeed, it is always possible to consider that the shear rate is expressed as $\dot{\gamma }=\frac{\dot{\gamma }}{1}$. If we now introduce the dimensionless shear rate $\dot{\widehat{\gamma }}=\frac{\dot{\gamma }}{{\dot{\gamma }}_1}$ where ${\dot{\gamma }}_1=1s^{-1}$, Equation (9) can be written:

$$\tau ={\tau }_K{\dot{\widehat{\gamma }}}^n,$$

(10)

where ${\tau }_K$ has the dimension of stress. It has exactly the same numerical value as K but its physical dimension is now independent of the consistency index, therefore explaining why these values can be compared with each other for different fluids.

The Herschel–Bulkley model combines both the Bingham plastic and the power law generalizations of the Newtonian rheological behavior [10]:

$$\tau ={\tau }_{\gamma }+{\tau }_K{\dot{\widehat{\gamma }}}^n.$$

(11)

It reduces to the Bingham plastic model when $n=1$ and to the power law model when ${\tau }_{\gamma }=0$. The Herschel–Bulkley model is also referred to as the yield power law (YPL) model.

Robertson and Stiff (1976) [11] proposed a model that has a smoother curvature than the YPL model:

$$\tau ={\tau }_A{\left(C+\dot{\widehat{\gamma }}\right)}^B,$$

(12)

where ${\tau }_A$ has the dimension of stress, and $B$, $C$ are dimensionless model parameters.

In medicine, biology and the food industry, the Heinz–Casson model [12] is often used:

$$\tau ={\left({{\tau }_{\gamma }}^n+{\left({\tau }_K\dot{\widehat{\gamma }}\right)}^n\right)}^{\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$n$}\right.},$$

(13)

as well as the Carreau model [13]:

$${\mu }_{eff}=\frac{\tau }{\dot{\gamma }}={\mu }_{\infty }+\left({\mu }_0-{\mu }_{\infty }\right){(1+{\left(\delta \dot{\gamma }\right)}^2)}^{\frac{n-1}{2}},$$

(14)

where ${\mu }_{eff}$ is the effective viscosity (dimension [ML⁻¹T⁻¹] (Pa.s)), ${\mu }_0$ is a Newtonian viscosity when the shear rate tends to zero, ${\mu }_{\infty }$ is a Newtonian viscosity when the shear rate tends to infinity, $\delta$ is a relaxation time (dimension [T] (s)) and $n$ is a power index.

Finally, Quemada developed a model [14] purely based on the physical principles established for the effects of deformable particles in suspension in a Newtonian fluid, which correspond rather well to typical drilling fluids where the viscosity behavior is caused by the level of entanglement of polymers or the assemblage of structures made of clay platelets:

$${\mu }_{eff}={\mu }_{\infty }{\left(\frac{1+{\left(\frac{\dot{\gamma }}{{\dot{\gamma }}_c}\right)}^n}{\sqrt{\frac{{\mu }_{\infty }}{{\mu }_0}}+{\left(\frac{\dot{\gamma }}{{\dot{\gamma }}_c}\right)}^n}\right)}^2,$$

(15)

where ${\dot{\gamma }}_c$ is a characteristic shear rate (dimension [T⁻¹] (s⁻¹)).

Each of these models can potentially describe the rheological behavior of drilling fluids. Yet, statistically, it can be shown that in a majority of cases, drilling fluids are better described by the Heinz–Casson, Quemada, Herschel–Bulkley, Robertson–Stiff and power law models in decreasing order, than by the Bingham plastic rheological behavior [15].

3.2. Pressure Losses in Laminar Flow

The viscous pressure gradient of the laminar flow of Newtonian fluid in a tube is well known. For Bingham plastic fluid, one can use the Swamee–Aggarwal [16] equation of the Darcy–Weisbach friction factor. As the Collins–Graves model does not have a yield stress, we can use the Weissenberg–Rabinowitsch–Mooney–Schofield method to estimate its associated laminar viscous pressure gradient. A solution is provided in Appendix A.

The laminar viscous pressure losses of a power law fluid in a pipe are also well-known [17]. The pressure loss of a Herschel–Bulkley fluid has been described by Kelessidis et al. (2006) [1], while the pressure loss of a Robertson–Stiff fluid is described in Beirute and Flumertfelt’s paper (1977) [18].

The viscous pressure loss gradient in laminar flow of a fluid following the Carreau rheological behavior can be efficiently calculated using the method described by Sochi [13], which is based on the Weissenberg–Rabinowitsch–Mooney–Schofield method. The same method can be used to estimate the laminar pressure gradient in a pipe of a fluid well-described by the Quemada rheological behavior. The derivation is given in Appendix B.

The pressure losses in laminar flow of a Casson fluid in a pipe have been described by Lee et al. (2011) [19]. Unfortunately, to the best of our knowledge, an analytic expression of the laminar pressure loss of a Heinz–Casson fluid in a pipe is not yet known. However, it is possible to perform a numerical estimation. The magnitude of the gradient of the fluid velocity at a radial position is given by the shear rate:

$$-\frac{dv}{dr}=\dot{\gamma }.$$

(16)

As we know the shear rate at the wall, since the shear stress at the wall is imposed by the pressure gradient, it is possible to integrate the fluid velocity in the tube from the wall, also considering that the velocity at the wall is zero as a consequence of the no slip at the wall hypothesis:

$$v_{j+1}=v_j-{\dot{\gamma }}_j\Delta r,$$

(17)

where $v_j$ is the fluid velocity at the radial position $r_j$, $v_{j+1}$ is the fluid velocity at the radial position $r_{j+1}=r_j-\Delta r$, $\Delta r$ is a radial step, ${\dot{\gamma }}_j$ is the shear rate at the radial position $r_j$, and considering that at $r_0=\frac{d_w}{2}$, $v_0=0$. Since ${\tau }_j=\frac{r_j}{2}\frac{dp}{ds}$, we can calculate ${\dot{\gamma }}_j$:

$${\dot{\gamma }}_j=\frac{1}{{\tau }_K}{\left({\left(\frac{r_j}{2}\frac{dp}{ds}\right)}^n-{{\tau }_{\gamma }}^n\right)}^{\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$n$}\right.}.$$

(18)

If ${\tau }_j$ gets smaller than ${\tau }_{\gamma }$, then the last value of the fluid velocity is used, here noted $v_{\gamma }$. The flowrate can then be estimated:

$$Q\approx {{\displaystyle \sum }}_{j=1}^n\pi v_j(r_j^2-r_{j+1}^2).$$

(19)

It is then simple to apply a Newton–Raphson method to find the value of $\frac{dp}{ds}$ that gives the desired volumetric flowrate. We have found that $n=100$ gives approximations of the pressure gradient with an accuracy better than 0.1%.

3.3. Pressure Gradient in a Tube in Transitional and Turbulent Flow

The Rabinowitsch–Mooney equation for a pipe is [20]:

$${\dot{\gamma }}_w={\dot{\gamma }}_{wN}\left(\frac{3}{4}+\frac{1}{4}\frac{d\ln \left(\frac{8\overline{v}}{d_w}\right)}{d\ln {\tau }_w}\right)=f_{RM}{\dot{\gamma }}_{wN},$$

(20)

where ${\dot{\gamma }}_w$ is the actual shear rate at the wall, ${\dot{\gamma }}_{wN}=\frac{8\overline{v}}{d_i}$ is the Newtonian shear rate at the wall, ${\tau }_w$ is the actual stress at the wall, $\overline{v}$ is the bulk fluid velocity and $f_{RM}$ is the Rabinowitsch–Mooney correction factor of the true shear rate at the wall compared to the Newtonian one. Note that for power–law and Herschel–Bulkley rheological behaviors, the Rabinowitsch–Mooney correction factor for pipe flow is:

$$f_{RM}=\frac{3}{4}+\frac{1}{4}n.$$

(21)

The Dodge and Metzner relation [21] can be used to estimate the Fanning friction factor ($f_F$) in a pipe for power law rheological behavior:

$$\frac{1}{\sqrt{f_F}}=\frac{4}{n^{0.75}}{\log }_{10}\left(R_e{f}_F^{1-\raisebox{1ex}{$n$}\!\left/ \!\raisebox{-1ex}{$2$}\right.}\right)-\frac{0.395}{n^{1.2}}.$$

(22)

Following the method described by Founargiotakis et al. [22], the actual shear stress at the wall can be approximated by a power law rheological behavior such that it follows locally the rheological behavior of the fluid around the given shear stress:

$${\tau }_w=K^{\prime }{\left({\dot{\gamma }}_{wPL}\right)}^{n^{\prime }},$$

(23)

where $K^{\prime }$ and $n^{\prime }$ are, respectively, the locally fitted consistency and flow behavior indices of the local power law rheological behavior and ${\dot{\gamma }}_{wPL}$ is the equivalent power law wall shear rate. From Equations (20) and (21), the expression of ${\dot{\gamma }}_{wPL}$ is:

$${\dot{\gamma }}_{wPL}=\left(\frac{3}{4}+\frac{1}{4}{n}^{\prime }\right){\dot{\gamma }}_{wN}.$$

(24)

The local power law flow behavior index can be calculated as follows:

$$n^{\prime }=\frac{d\ln {\tau }_w}{d\ln {\dot{\gamma }}_{wN}},$$

(25)

and the local power law consistency index is simply:

$$K^{\prime }=\frac{{\tau }_w}{{\left({\dot{\gamma }}_{wPL}\right)}^{n^{\prime }}}.$$

(26)

It is then possible to define a local power law Reynolds number ($R_e^{\prime }$) as:

$$R_e^{\prime }=\frac{\rho {\overline{u}}^{2-n^{\prime }}{d}_w^{n^{\prime }}}{K^{\prime }{\left(\frac{2}{n^{\prime }}+6\right)}^{n^{\prime }}},$$

(27)

where $\rho$ is the fluid density. The local power law Reynolds number limits for laminar ($R_{e_l}^{\prime }$) and turbulent ($R_{e_t}^{\prime }$) flow are [21]:

$$R_{e_l}^{\prime }=3250-1150n^{\prime },$$

(28)

$$R_{e_t}^{\prime }=4150-1150n^{\prime }.$$

(29)

The Darcy–Weisbach friction factor at the laminar limit ($f_{lDW}$) is:

$$f_{lDW}=\frac{64}{R_{e_L}^{\prime }},$$

(30)

while the Fanning friction factor at the limit of turbulent flow ($f_{tF}$) is estimated using Equation (22) with $R_{e_t}^{\prime }$. Note that the Darcy–Weisbach friction ($f_{DW}$) is simply related to the Fanning friction ($f_F$) by:

$$f_{DW}=4f_F.$$

(31)

For transitional flow, the Darcy–Weisbach friction ($f_{trDW}$) is linearly interpolated between the laminar and turbulent limits as a function of the local power law Reynolds number:

$$f_{trDW}=f_{lDW}+\frac{R_e^{\prime }-R_{e_l}^{\prime }}{R_{e_t}^{\prime }-R_{e_l}^{\prime }}\left(f_{tDW}-f_{lDW}\right).$$

(32)

The pressure gradient for transitional and turbulent flow is then calculated using the Darcy–Weisbach definition of the friction factor:

$$\frac{dp}{ds}=f_{DW}\frac{\rho {\overline{u}}^2}{2d_w}.$$

(33)

3.4. Comparison of Estimated Pressure Losses for Different Drilling Fluids

In this section we will compare the influence of the matched rheological behavior with expected pressure losses in a circular pipe. All the examples are based on measurements made on real drilling fluids with a scientific rheometer and utilizing rotational speeds spread along a logarithmic scale in order to obtain a proper coverage at low shear rates.

The example of Figure 3 corresponds to a KCl/polymer water-based mud of density 1250 kg/m³ at 80 °C (composition in % by wt: fresh water 65.35, KCl 9.60, Soda ash 0.08, DuoTec NS 0.32, Trol FL 0.80, Glydril MC 2.40, barite 19.45). We can see that the best match is obtained with a power law rheological behavior, or equivalently, a Robertson–Stiff, Herschel–Bulkley or Heinz–Casson models as they can all degenerate to a power law formulation.

When applying the various rheological behaviors to pressure loss calculations in a 5 in pipe for flowrates ranging between 0 and 5000 L/min, we can see that the expected pressure gradients are rather different for the Collins–Graves, Quemada and Heinz–Casson models than for the power law rheological behavior (see Figure 4). For the Quemada and Heinz–Casson models, the difference in pressure gradient (relative to the power law model) reaches 3% in laminar flow and 29% in turbulent flow, even though the flow curve error with the actual rheometer measurements does not exceed −2.5% on the under-evaluation side and 11% on the over-evaluation side.

Another KCl/polymer fluid (composition in % by wt: fresh water 50.92, KCl 8.00, soda ash 0.07, DuoTec NS 0.27, Trol FL 0.67, Glydril MC 2.00, barite 38.07) with a density of 1500 kg/m³ at 80 °C, is best described by the Herschel–Bulkley or the Robertson–Stiff model (see Figure 5). The ${\chi }^2={\displaystyle \sum }_{i=1}^n{\left({\widehat{\tau }}_i-{\tilde{\tau }}_i\right)}^2$ for the Herschel–Bulkley model is 0.013, while for the Robertson–Stiff model, it is 0.014. Note that here ${\widehat{\tau }}_i$ is the measured shear stress and ${\tilde{\tau }}_i$ is the modelled shear stress. The Robertson–Stiff model fits best at low shear rates, while the Herschel–Bulkley model is better at high shear rates. Other models such as Heinz–Casson, Quemada and Carreau are quite close as well (max percentage difference 7% compared to the YPL model). However, the power law and Collins–Graves models do not match the rheometer measurements so well. Yet, when simulating pumping the fluid in a 5 in pipe, the estimated pressure gradients are rather dispersed, with differences in the range of 2% in laminar flow, but exceeding 35% toward the end of the analyzed flowrates for the turbulent flow regime, each of the models giving a relatively different estimation of the pressure gradient as a function of the flowrate (see Figure 6).

In practice, depending on the fluid, any of the above-listed rheological behaviors may provide the best match with rheological measurements. For instance, a KCl/polymer with a density of 1750 kg/m³ at 50 °C (composition in % by wt: fresh water 39.18, KCl 6.86, Soda ash 0.06, DuoTec NS 0.23, Trol FL 0.57, Glydril MC 1.71, barite 51.39) may be best modeled by a Quemada model, while the Collins–Graves model may provide the best fit for a bentonite-based fluid as it can be seen in Figure 7 (composition in % by wt: fresh water 94.23, bentonite 5.49, caustic soda 0.27, soda ash 0.01).

From the above described examples, it seems that the predicted pressure losses, especially in turbulent flow, can differ significantly as a function of the chosen rheological model. Therefore, it should be a good practice to analyze which model is best suited for a fluid rather than choosing one model that should fit all possible fluids. This being said, it remains to be verified that finding the best rheological model that matches rheometer values helps improving the pressure gradient predictions.

3.5. Verification in a Laboratory Flow-Loop

To verify whether the choice of the best fitting rheological behavior based on rheometer measurements is also discernable when comparing actual pressure gradients in a pipe and estimated values, as described in Section 3.2 and Section 3.3, we recourse to measurements made in a laboratory flow-loop. A positive displacement pump circulates fluid stored in a thermo-regulated tank with internal circulation. First the fluid flows into a pipe with an internal diameter of 25 mm and then into a glass tube with an internal diameter of 15.5 mm. Along the glass tubes there are three differential pressure sensors. The distances between the ports of the differential pressure sensors are, respectively, 0.209 m, 0.212 m and 0.206 m, following the direction of the flow. The first sensor is placed 1.955 m away from the change in diameter from 25 mm to 15.5 mm. The second sensor is positioned 1.415 m apart from the first one and the last sensor is separated from the previous by 1.46 m (see Figure 8). There is also a pressure reducing valve placed 1.5 m upstream from the change in diameter. The purpose of this valve is to ensure that the fluid is fully sheared before it continues its journey inside the flow-loop.

All the experiments were conducted at a controlled temperature of 25 °C. The temperature variations during the experiments have not exceeded ±0.5 °C. The first measurements have been made with an aqueous solution of glycerin (50 %vol) having a density of 1115 kg/m³ at 25 °C. This is a Newtonian fluid with a viscosity of 0.00445 Pa.s (see Figure 9a).

Figure 9b shows the measured pressure gradients taken with the three differential pressure sensors. The measurements from the first differential pressure sensor have been filtered out for Reynolds number above 4000 as entrance effects influenced the readings. This corresponds to a turbulent flow regime for which entrance effects can perturb measurements at a larger distance than in laminar flow [23]. The calculated pressure gradient is shown as the solid brown curve. It is in good agreement for laminar and turbulent flow, but does not match well in the transitional flow regime. These results are essentially provided for verification purposes.

The second experiment makes use of an aqueous solution of Carbopol with a low molecular weight at the concentration of 0.07%. The acidity of the solution is neutralized by adding sodium hydroxide in a concentration of 0.009%. The fluid is shear-thinning and well-described by a power law rheological behavior (${\tau }_K=0.019 \mathrm{Pa}.\mathrm{s}, n=0.848$), seconded by the Quemada model (${\mu }_{\infty }=0.0053 \mathrm{P}\mathrm{a}.\mathrm{s}, {\mu }_0=0.1 \mathrm{P}\mathrm{a}.\mathrm{s},{\dot{\gamma }}_c=14 {\mathrm{s}}^{-1}, n=0.35$) and the Heinz–Casson rheological behavior (${\tau }_{\gamma }=0.08 \mathrm{P}\mathrm{a}, {\tau }_K=0.046 \mathrm{P}\mathrm{a}, n=0.50$) (see Figure 10a). In Figure 10b, we can see that for this particular fluid, the power law rheological behavior provides more accurate pressure gradients for the whole range of measurements than the estimations made with the Quemada and Heinz–Casson models.

However, for the laminar flow regime, pressure loss gradients estimated with different rheological models do not differ much. Figure 11a shows the measured rheological behavior of an aqueous solution of polyanionic cellulose (concentration 0.2%) measured at 25 °C with a scientific rheometer together with curve-fitted models for power law (${\tau }_K=0.070 \mathrm{P}\mathrm{a}.\mathrm{s}, n=0.728$), Heinz–Casson (${\tau }_{\gamma }=0.19 \mathrm{P}\mathrm{a}, {\tau }_K=0.069 \mathrm{P}\mathrm{a}, n=0.324$) and Quemada (${\mu }_{\infty }=0.00323 \mathrm{P}\mathrm{a}.\mathrm{s}, {\mu }_0=0.100 \mathrm{P}\mathrm{a}.\mathrm{s},{\dot{\gamma }}_c=1845 {\mathrm{s}}^{-1}, n=0.35$) models. As it can be seen in Figure 11b, the estimated pressure gradients are almost indistinguishable from each other. Note that past a volumetric flowrate of 21 l/min, the flow regime is transitional.

However, it should be noted that transitional flow is observable on the differential pressure measurements, by the appearance of spikes of pressures. Just past the upper limit of laminar flow, these spikes are only occasional and of low amplitudes, but the more the Reynolds number increases within the transitional flow regime, the more frequent the spikes are and the larger are their amplitudes (see Figure 12). This effect is not reflected through the use of an interpolation of the pressure gradient between the laminar and turbulent pressure losses as described by Equation (32), as the result is only deterministic while in reality the phenomenon is stochastic in nature.

When fitting the two above examples with a Herschel-Bulkley or a Robertson-Stiff rheological behavior, the fitted value of the yield stress turns out to be zero. If we utilize a zero yield stress in the Heinz-Casson rheological behavior, the model reduces to the one of a Newtonian fluid: $\tau ={\left(0^n+{\left({\tau }_K\dot{\widehat{\gamma }}\right)}^n\right)}^{\raisebox{1ex}{$1$}\!\left/ \!\raisebox{-1ex}{$n$}\right.}={\tau }_K\dot{\widehat{\gamma }}$ and therefore is not capable to reproduce the shear-thinning tendency that is measured with the rheometer. Therefore, the fitted parameters of the Heinz-Casson model always utilize a non-zero yield stress as soon as there is a shear-thinning tendency, even though, other models like Herschel-Bulkley or Robertson-Stiff would indicate that the yield stress is zero.

Figure 13a shows the measured rheogram of an aqueous solution of Carbopol with high molecular weight. This time the fluid has a yield stress. The Herschel-Bulkley parameters are: ${\tau }_{\gamma }=1.21 \mathrm{P}\mathrm{a}, {\tau }_K=0.266 \mathrm{P}\mathrm{a}, n=0.642$, the parameters for the Robertson-Stiff model are: ${\tau }_A=380 \mathrm{P}\mathrm{a}, B=0.597, C=10.458$, those of the Heinz-Casson model are: ${\tau }_{\gamma }=0.24 \mathrm{P}\mathrm{a}, {\tau }_K=0.488 \mathrm{P}\mathrm{a}, n=0.271$ and the parameters of the Quemada model are: ${\mu }_{\infty }=0.00604 \mathrm{P}\mathrm{a}.\mathrm{s}, {\mu }_0=\mathrm{10,000} \mathrm{P}\mathrm{a}.\mathrm{s},{\dot{\gamma }}_c=970 {\mathrm{s}}^{-1}, n=0.364$. The comparison between the estimated pressure gradients and the measured ones is shown in Figure 13b. There are no significant differences between the predicted pressure gradients made with the different models. There is a good agreement between the measured pressure gradients and the estimated ones made with the rheological models calibrated with the rheometer measurements, until approximatively 22 L/min. However, for the largest flowrates, there is a discrepancy building up, as the measurements tend to indicate that the fluid becomes less shear-thinning than what is expected from the rheometer measurements.

All the previous fluids are non-thixotropic. This can be verified by confirming the time independence of the shear stress after step-changing the rheometer speed. Figure 14 shows a series of measurements made with a scientific rheometer where the Newtonian shear rate is changed from 1020, 512, 341, 170, 102, 51, 10.2, 5.1, 10.2, 51, 102, 341, 512, 1020, consecutively. It is noticeable that the polyanionic cellulose high viscosity (PAC-HV) 0.2% Carbopol with low molecular weight and the second aqueous solution of Carbopol with high molecular weight respond immediately to the change in shear rate.

Figure 15 shows that with an aqueous solution of KCl and xanthan gum, the shear stress stabilizes to a constant value only after several tens of seconds after the shear rate has changed. The KCl/polymer fluid is therefore thixotropic.

Figure 16a shows the steady state rheogram measured with a KCl/xanthan gum fluid while Figure 16b shows that the three differential pressure sensors measure different pressure gradients until a flowrate of 25 L/min, in contradiction with what has been observed when utilizing non-thixotropic fluids.

As it seems that thixotropic effects influence the pressure gradients at different positions along a pipe after a change in diameter, we will analyze, in the next section, how to incorporate thixotropy in pressure loss calculations for flow in circular pipes.

4. Laminar Flow Pressure Losses of a Thixotropic Fluid in a Hydraulic Circuit with Change in Diameters

The thixotropic response of a fluid is described by a dimensionless variable called the structure parameter ($\lambda$). The structure parameter corresponds to the degree of entanglement of floccules in the fluid, such as that generated by the buildup of structures by polymers or clay platelets. Mewis and Wagner (2009) [24] propose a general description of the time evolution of the structure parameter through the following ordinary differential equation:

$$\frac{d\lambda }{dt}=k_1{\dot{\gamma }}^a{\left(1-\lambda \right)}^b-k_2{\dot{\gamma }}^c{\lambda }^d,$$

(34)

where $k_1$ and $k_2$ are, respectively, the coefficients that characterize the breakdown and the buildup of floccules, and $a,b,c$ and $d$ are either explicitly specified by the model or fitted to a series of observations. During the derivation of his model, Quemada [14] utilizes $a=0$, $b=1$, $c=1$ and $d=1$ and we will follow the same assumption. The structure parameter is then used to modify the flow curve. Here, we will utilize a simple first order function of the structure parameter:

$$\{\begin{array}{c}{\tau }_{\dot{\gamma },\lambda }=\left(1+\alpha \frac{\lambda -{\lambda }_{\dot{\gamma },\infty }}{{\lambda }_{\dot{\gamma },\infty }}\right){\tau }_{\dot{\gamma },\infty }\\ \frac{d\lambda }{dt}=k_1\left(1-\lambda \right)-k_2\dot{\gamma }\lambda \end{array},$$

(35)

where $\alpha$ is a scaling factor, ${\lambda }_{\dot{\gamma },\infty }$ is the value of the structure parameter for the shear rate $\dot{\gamma }$ in steady state conditions, i.e., when $t\to \infty$, and ${\tau }_{\dot{\gamma },\infty }$ is the best match rheological behavior in steady state conditions, e.g., power law, Collins–Graves, Herschel–Bulkley, Robertson–Stiff, Heinz–Casson, Carreau, Quemada. The integration of the ordinary differential equation of the structure parameter at constant shear rate gives:

$$\lambda =\beta e^{-\left(k_1+k_2\dot{\gamma }\right)t_{step}}+\frac{k_1}{k_1+k_2\dot{\gamma }},$$

(36)

where $\beta$ is an integration constant and $t_{step}$ is the time since the shear rate changed. If we assume a continuity of the structure parameter when the shear rate changes and if we consider that the origin of time is set at the moment when the shear rate is modified, then:

$$t_{step}=0, {\lambda }_0=\frac{k_1}{k_1+k_2\dot{\gamma }}+\beta \iff \beta ={\lambda }_0-\frac{k_1}{k_1+k_2\dot{\gamma }},$$

(37)

denoting ${\lambda }_0$ the value of the structure parameter when the shear rate changed. The value of ${\lambda }_{\dot{\gamma },\infty }$ can then be calculated:

$${\lambda }_{\dot{\gamma },\infty }=\underset{t_{step}\to \infty }{\lim }\lambda =\frac{k_1}{k_1+k_2\dot{\gamma }}.$$

(38)

When the flowrate is constant, there is no acceleration in the portions of pipes with a constant diameter and, therefore, Equation (1) and consequently Equation (3) still yield, regardless of the thixotropic behavior. The hydraulic circuit is now discretized in short sections of length $\Delta s$ (dimension [L] (m)) and each of these sections is indexed by the subscript $i$. The cross-section is discretized in $n$ small shells of radial thickness $\Delta r_i=\frac{d_{w_i}}{2n}$ (dimension [L] (m)), where $d_{w_i}$ is the internal diameter at along pipe position $i$. The radial position is indexed with the subscript $j$, starting at zero at the outer position in the cross-section. In laminar flow, we assume that the fluid contained in the shell $i-1,j$ moves to the shell $i,j$ after an elapsed time $\frac{v_{i-1,j}}{\Delta s}$, where $v_{i-1,j}$ is the fluid velocity at along hole and radial position $i-1,j$ (see Figure 17). It may be that the diameter of the pipe changes from along tube position $i-1$ to $i$ from a diameter $d_{w_{i-1}}$ to a diameter $d_{w_i}$. Assuming that the fluid is incompressible within the working pressure, the mass conservation equation can simply be expressed as a volume conservation:

$$\begin{array}{c}{v}_{i-1,j}\left({\left(\frac{d_{w_{i-1}}}{2}-j\frac{d_{w_{i-1}}}{2n}\right)}^2-{\left(\frac{d_{w_{i-1}}}{2}-\left(j+1\right)\frac{d_{w_{i-1}}}{2n}\right)}^2\right)=v_{i,j}({\left(\frac{d_{w_i}}{2}-j\frac{d_{w_i}}{2n}\right)}^2-\\ {\left(\frac{d_{w_i}}{2}-\left(j+1\right)\frac{d_{w_i}}{2n}\right)}^2),\end{array}$$

(39)

and therefore, we can express the fluid velocity at the entrance of the shell with axial index $i$:

$$v_{i,j}=\frac{v_{i-1}\left({\left(\frac{d_{w_{i-1}}}{2}-j\frac{d_{w_{i-1}}}{2n}\right)}^2-{\left(\frac{d_{w_{i-1}}}{2}-\left(j+1\right)\frac{d_{w_{i-1}}}{2n}\right)}^2\right)}{\left({\left(\frac{d_{w_i}}{2}-j\frac{d_{w_i}}{2n}\right)}^2-{\left(\frac{d_{w_i}}{2}-\left(j+1\right)\frac{d_{w_i}}{2n}\right)}^2\right)}.$$

(40)

It is then possible to perform a numerical integration along the hydraulic circuit, considering that the initial conditions of the structure parameter at each radial position are imposed by the fact that the fluid has been flowing through a long section of pipe with a constant diameter and therefore the structure parameter is equal to ${\lambda }_{\dot{\gamma },\infty }$ anywhere in a cross-section upstream from the change in diameter. The same methodology is followed as for solving the pressure losses with a Heinz–Casson fluid, except that we need to keep track of the time evolution of the structure parameter for each radial position and at each axial position.

In the laboratory flow-loop, when pumping at a constant volumetric flowrate of 6 L/min in the example flow-loop depicted in Figure 8, with a thixotropic fluid where $k_1=0.003$, $k_2=0.001$, $\alpha =0.05$, the structure parameter is impacted by the change in pipe diameter from 25 to 15.5 mm (see Figure 18). It is mostly the structure parameter closest to the pipe wall which is affected by the change in geometry. In the central region of the pipe, there is an unsheared region corresponding to the yield stress of the fluid. This part is not much influenced after the transition from the larger to the smaller diameter.

The axial evolution of the structure parameter impacts in turn the relationship between the shear stress and shear rate. As the shear stress needs to be a linear function of the radial position (to the limit of yield stress), this impacts the acceptable shear rate and consequently how the structure parameter evolves in the axial direction. Figure 19 shows the estimated values of the shear rate in both the axial and radial directions. We can see that the change in diameter continues to impact the shear rate over a long distance after the change in diameter.

As a consequence, the pressure gradient after the change in diameter is larger than if it were estimated with a non-thixotropic fluid having a similar rheological behavior. The pressure gradient decays to the one corresponding to steady state conditions at infinite distance to any change in geometrical dimensions (see Figure 20).

The difference of the pressure gradient estimated with a thixotropic fluid compared to a non-thixotropic fluid with equivalent rheological behavior in steady state conditions increases with the volumetric flowrate (see Figure 21). The difference is of course largest the closest to the change in diameter, but it is estimated to still have effects more than 5.5 m after the change in diameter in the context of the laboratory flow-loop. With that particular fluid and with the described flow-loop dimensions, it is necessary to position a differential pressure sensor at least 7.7 m downstream from the change in pipe diameter to be within 1% of the pressure gradients for which thixotropic effects would be negligible.

5. Discussion

The estimation of pressure losses for drilling operations can basically be split into the estimation of pressure losses inside the string and in the annulus. The paper focuses on the estimation of pressure losses in circular pipes and, therefore, is applicable to pressure losses inside the string. Such estimations are important when determining the maximum flowrate possible to use when being limited by the maximum pressure rating of the surface installation pumping system. Such pressure gradients are also important when estimating the dynamic response of the drilling system to pump accelerations and more generally when solving the Navier–Stokes equations for the complete hydraulic circuit.

In practice, it is not uncommon that the estimated pump pressures made with real-time hydraulic calculations are off by several tens of bars. There may be many reasons for such inaccuracies, including, but not limited to, the lack of precision of rheometer measurements made at the rig site, the accuracy of pressure and temperature extrapolation models for the fluid rheological behavior, uncertainty related to the actual pipe internal diameter, unprecise modelling of the effect of change in diameters at the level of tool-joints. However, as presented above, simply the choice of the rheological model may substantially influence the estimated pressure losses in the turbulent flow regime. Indeed, most of the time, the flow regime is turbulent inside the pipes during normal drilling operations. Yet, there is no special information available that may help decide which rheological model is the most pertinent to use for the turbulent flow regime.

This problem could be resolved if inline pipe rheometers performed measurements in the turbulent flow regime. Such pipe rheometers could indicate which rheological behavior is best suited to estimate pressure losses in pipes when subjected to turbulent flow. At the same time, they could provide the parameters of the best-suited rheological model. The calibration of these parameters requires one to apply low flowrates in order to accurately estimate the yield stress, if any, of the drilling fluid. Yet, we have seen in Section 4 the importance of the change in diameters on the estimated pressure gradients in a pipe rheometer when utilizing thixotropic fluids. The estimated parameters of the rheological model may become biased by the impact of structuration and de-structuration of floccules after passing a restriction or an enlargement. It is only past long distances that these effects become negligible and therefore the size of a pipe rheometer could turn out to be unreasonably large to eliminate those effects. Another solution would be to position several differential pressure sensors, as shown in Figure 8, in order to collect information about thixotropic effects and thereafter calibrate the thixotropic model described by Equation (35).

As a continuation of this work, we will, in the near future, investigate the impact of time dependent changes on the boundary conditions, for instance oscillations of the volumetric flowrate, when applied to thixotropic fluids. Figure 22 illustrates initial measurements that have been made in the laboratory flow-loop. It seems that the time dependence of the rheological properties of the thixotropic fluid influences the rapidity by which the fluid can respond to quick changes in flow conditions. This is an important step to study before investigating the impact of thixotropy on pressure losses in an annulus under the influence of a moving inner pipe, both axially, radially and laterally.

6. Conclusions

We have derived semi-analytical solutions for the estimation of pressure gradients in cylindrical pipes with the Collins–Graves and Quemada rheological behaviors in laminar flow conditions. A method to numerically calculate the pressure gradients with the Heinz–Casson rheological behavior has also been described. It has served as a basis for estimating pressure gradients with thixotropic fluids. A general method has been explained to numerically calculate pressure gradients in turbulent flow.

We have seen that different rheological behaviors, fitting accurately with high precision rheometer measurements, give similar pressure gradients in the laminar flow regime. However, the differences may be substantial when estimating the pressure gradients in turbulent flow.

We have also seen that for thixotropic fluids, the passage through a change in diameter may influence the pressure gradients over substantial distances. We have exposed a method to estimate the pressure gradients in steady state conditions for thixotropic fluids circulating in a hydraulic circuit containing diameter variations.

The results have been confirmed by measurements made with a laboratory scale flow-loop.