The rheological behavior of energized fluids and foams with application to hydraulic fracturing: Review

1. Introduction

In the last decade, resources such as shale gas and tight oil have become extremely important thanks to the successful development of hydraulic fracturing (Smith and Montgomery, 2015). The resulting increase in the use of natural gas has allowed the United States economy to prosper, even while reducing its carbon footprint. Despite these positive aspects, hydraulic fracturing is controversial and widely threatened by regulations, to the point where it has been banned in several areas. The reasons for this include the consumption of large amounts of water and the risks inherent in its disposal. The productivity of hydraulically fractured wells also declines very quickly, and is usually much less than would be expected considering the resources in the reservoir. Thus, there is a strong motivation to improve the productivity and lifetime of fractured wells, while reducing the quantity of water and other chemicals consumed in the fracturing process.

Foam fracturing fluids are promising in both respects, as a significant fraction of water is replaced by nitrogen or carbon dioxide. The use of foam fracturing fluids is well-established in the small percentage of important reservoirs that are under-pressured, where the quick and easy clean-up offered by foams is sufficient to justify the effort, specialized equipment and expertise required. Their wider use is hampered by logistical challenges, mainly because modern hydraulic fracturing is increasingly applied in long horizontal wells with 20–40 stages which require delivery of several hundred truckloads of liquid nitrogenover a few days, or even higher volumes in the case of carbon dioxide (because of the higher density of carbon dioxide under downhole conditions). However, these challenges are manageable and foams can provide productivity benefits beyond under-pressured reservoirs, mainly due to improved proppant transportand reduced water damage, see e.g., (Palmer and Sito, 2013, McAndrew et al., 2014a). In this review, we will discuss the basis for these benefits, and particularly how recent scientific insights support their implementation.

In hydraulic fracturing, pressurized fluids are used to fracture impermeable rockand transport proppant into the fracture to prevent fracture closure after the pressure release. Proppant placement within the fracture directly impacts productivity, because it controls both short-term and long-term conductivities of the fractured well. To prevent fracture closure and enhance the conductivity of wells, different fracturing fluids and proppant types have been used, see review papers by Al-Muntasheri (2014) and Liang et al. (2015), respectively. The main types of fluids are: (i) low viscosity fluids (e.g. slickwater, which is water with a drag-reducing additive at a few percent by weight) that are efficient in creating long fractures in the flow direction, but not in proppant transport due to a high rate of proppant segregation (Palisch et al., 2010), (ii) gels and polymer-based fluids that lead to fractures that propagate less in the flow direction, more in the direction perpendicular to the flow, and transport proppant well, but may damage the fracture surface (Gomaa et al., 2014), (iii) energized fluids (bubbly liquids) and foams that have a high effective viscosity and therefore behave more like gels, but avoid damage to the fracture surface and reduce water consumption (Palmer and Sito, 2013, Wanniarachchi et al., 2015), and (iv) a mixture of the aforementioned fluids using methods like hybrid (Sharma et al., 2004), reverse hybrid (Liu et al., 2007), multi-stage (Manchanda and Sharma, 2013), alternate-slug (Malhotra et al., 2014) fracturing.

In under-pressurized geologic plays, energized fluids (see definition in section 2) and foams serve better than other fluids in all steps of the treatment, starting from drilling (Loureno et al., 2004, Kelessidis et al., 2006, Edrisi et al., 2014), fracturing (Palmer and Sito, 2013), proppant transport (Valko and Economides, 1997, Deshpande et al., 2013, Gomaa et al., 2015b, Zhang, 2014), and finally deep well cleaning (Gandossi, 2013). In terms of geography, the use of energized fluids and foams is most prevalent in Canada and in severely water-stressed regions of the U.S. (mainly in New Mexico). Multiple studies, see e.g., Burke et al., 2011, Reynolds et al., 2015, have been published discussing the performance of foam fracturing fluids in the Montney and Cardium (in Canada). Again, these reservoirs are largely under-pressured. A few studies, mainly from the 1980's- 1990's, have discussed foam fracturing applications in more highly-pressured reservoirs. For example, Harris et al. (1984). discussed results from the Red Fork formation, which is at a depth of 3900 m and highly over-pressured. A larger study of over 85 wells was carried out in strata above the Haynesville (Warnock et al., 1985), with depths to below 4300 m and pressures as high as 910 bar, some of the greatest depths ever recorded for foam fracturing.

Polymeric and polymer-free foams (or energized fluids) are viscoelastic materials whose mechanical and rheological responses to transmitted stresses fall between elastic solids and viscous fluids (Zhang et al., 2015, Gomaa et al., 2015a). The mechanical and rheological responses and flow patterns of such fluids are strongly affected by the stress conditions (e.g. applied shear rates), flow geometry (e.g. pipe or fractures size), time-scale (e.g. aging and drainage), microstructures (e.g. bubble size distribution and volume fraction, packing configuration, surfactant molecules, and other chemical additive), and state physical parameters (e.g. operating temperature and pressure).

During the initial development of foam fracturing (1970s–1990s), multiple studies (e.g (Harris et al., 1984, Princen and Kiss, 1989).) described their rheology from a more or less practical point-of-view and with particular application to fracturing. Since that time, fundamental understanding of the complex rheological behavior of foams has made considerable progress through many studies, covering aqueous foams (Herzhaft, 1999, Stevenson, 2012, Dollet and Raufaste, 2014, Fameau et al., 2015), non-aqueous and polymeric foams(Blazquez et al., 2014), polymer-free foams (Herzhaft, 1999, Gu and Mohanty, 2015), and particle-laden stabilized foams (Emrani and Nasr-El-Din, 2015). In this review we will point out relatively new fundamental concepts that are particularly relevant to the use of foams in fracturing, notably the following:

•
High bubble deformation in narrow fractures and near the fracture tip, which will reduce the apparent viscosity and facilitate fracture growth (Faroughi and Huber, 2015).
•
Osmotic pressure, or the tendency of foams to retain water driven by interfacial energies. Thanks to this effect, leakoff of water from foams is much lower than would be expected based on the volume fraction of water present. Lower leakoff from foams has been observed in the laboratory but its origin in osmotic pressure was not explained (Ribeiro and Sharma, 2012).
•
Elastic behavior of foams which can improve proppant transport without necessarily requiring an increase in apparent viscosity (Gomaa et al., 2014, Gomaa et al., 2015a, Gomaa et al., 2015b).

This review is organized as follows. In Section 2, we explain the differences between energized fluids and foams, and explain the effect of the dispersed phase volume fraction on the properties of each. In sections 3 Energized fluids, 4 Foams, we present deeper reviews of the properties of energized fluids and foams, respectively. Section 5 discusses how these concepts are applied in fracturing. Next, in section 6, we conclude with a discussion of the potential for broader application of foams in fracturing and highlight some areas where work is ongoing and more understanding is needed.

2. Regimes of liquid-gas mixture

A fluid-fluid mixture is a system where two or more fluids with different physical properties are mixed, but they are not combined chemically (i.e. they are immiscible). In these systems, one fluid is continuous, called the ambient phase, and the rest are suspended, called dispersed phase(s). Depending on the nature of the dispersed phase(s), the fluid is classified either as an emulsion (liquid) or a foam (gas). In the petroleum literature, an additional distinction is often made within the latter category between mixtures with low volume fraction of gas dispersed, called ”energized fluids”, while the label ”foam” is restricted to high gas volume fraction systems. The term ”energized fluid” is somewhat confusing, especially to those not in the petroleum industry. Its origin is in the use of such fluids to provide ready flowback from the well once pressure is released, even in an underpressured reservoir (Obo, 2015) - the fluid is said to have ”energy” because it flows back by itself, without requiring reservoir pressure or a pump. Of course, this benefit applies to foams also, so that strictly speaking ”energized fluid” should be a generic term that includes foams. However, in practice it is reserved for fluids where the gas volume fraction is too low to provide a sharp increase in apparent viscosity, similar to what is usually called a ”bubbly liquid” in the foam physics literature. We will follow the petroleum literature terminology here.

In a stress-free state, the dispersed phase typically forms spherical bubbles (energized fluids and foams) or droplets (emulsions) inside the continuous fluid. Hereafter, we refer to the dispersed phase as ”bubbles” in all cases. In our case of interest in this review, a gas (typically nitrogen, $N_{2}$ or carbon dioxide, $C O_{2}$ ) is dispersed inside an aqueous phase forming an energized fluid or a foam. Fig. 1 shows a schematic two-dimensional representation of a two-phase gas-liquid mixture and the classification based on the bubble volume fraction (i.e. dispersed phase volume fraction, ψ). This parameter is known as ”quality” in the petroleum industry. Different categories based on volume fraction can be identified at this point.

Within the energized fluid category, at a very low volume fraction known as dilute regime, $0 \leq ψ \leq 0.1$ , bubbles are few and far between in the ambient fluid (see Fig. 1). Collisions between bubbles are almost nonexistent, and the assumption of no hydrodynamic interactions between bubbles is valid (Einstein, 1911, Batchelor and Green, 1972, Brady, 1993, Lim et al., 2004). The effective rheological behavior of dilute energized fluids is a weak linear function of the bubble volume fraction following Einstein, 1911, Taylor, 1932, Mackenzie, 1950, Frankel & Acrivos (Frankel and Acrivos, 1970) when $ψ \to 0$ , and a weak non-linear function following the work of Saito, 1950, Mendoza, 2011, Faroughi & Huber (Faroughi and Huber, 2015) in the entire range of $0 \leq ψ \leq 0.1$ .

With the increase of dispersed phase volume fraction, the disturbance of the background flow streamlines due to presence of bubbles increases leading to a weak hydrodynamic interaction between bubbles (see Fig. 1 where these interactions are symbolized with dashed arrows). This regime is known as semi-dilute regime where $0.1 \leq ψ \leq 0.2$ (Choi and Schowalter, 1975, Goddard and Miller, 1967, Leal-Calderon et al., 2007a). When bubbles start to interact, mediated by the ambient fluid, the rate of viscous dissipation may increase depending on the shear condition and bubble deformation. These interactions thus invalidate the linear relationship between the relative viscosity and bubble volume fraction (Oldroyd, 1959, Rust and Manga, 2002, Lim et al., 2004).

For regimes with bubble volume fraction of $0.2 \leq ψ \leq 0.4$ , known as intermediate regimes, collisions between bubbles become possible, and strong hydrodynamic interactions develop between bubbles. The force balance in this regime is still dominated by the shear and surface tension stresses at the interface of mixture constituents. Many experimental studies show the presence of slight shear thinning behavior (reduction in relative viscosity with increasing shear rate) due to the shear-induced microstructural rearrangement. At a given shear condition, the effect of bubble volume fraction on the relative viscosity becomes more pronounced within this regime (Ducamp and Raj, 1989, Manga and Loewenberg, 2001).

Beyond 40% volume fraction of bubbles, the rheological and mechanical behavior of energized fluid becomes very complex (Eilers, 1943, Pal, 2000, Pal, 2001, Pal, 2003, Llewellin and Manga, 2005, Derkach, 2010, Faroughi and Huber, 2015). This degree of complexity arises due to high frequency of bubble-bubble interactions as well as the complex many-body hydrodynamic interactions occurring between bubbles. This regime ( $0.4 \leq ψ \leq ψ_{t}$ where $ψ_{t}$ is a threshold bubble volume fraction at which the maximum possible volume of the dispersed phase is made of spherical bubbles) is known as the concentrated regime, where bubbles are placed close and form common boundaries (i.e. made of films of continuous phase) together. For suspensions of solid spherical particles, $ψ_{t}$ is known as the amorphous maximally random jammed packing at which the transition between fluid-like and solid-like states occurs, and generally found to vary in $0.58 < ψ_{t} \leq 0.637$ for a monomodal particle size distribution (Scott and Kilgour, 1969, Song et al., 2008, Boyer and GuazzelliPouliquen, 2011). For any mixture, the value of the threshold packing volume fraction depends strongly on the state of stress, bubble size distribution and the packing configuration, i.e. the presence or absence of randomness, and thus in energized fluids with spherical monomodal particles, $ψ_{t}$ may vary between 0.52 (corresponding to cubic lattice packing) to 0.74 (corresponding to the face centered cubic or hexagonal close packing). We note that the value of $ψ_{t}$ increases as modality (e.g. from monomodal to bimodal) increases (Faroughi and Huber, 2014). The interactions between bubbles in this regime are both direct and indirect, which changes the rate by which momentum is exchanged between constituents of the mixture. Other phenomena like lubrication and the gradient of surface tension caused by uneven distribution of surfactant concentration also affect the bulk properties of concentrated energized fluids. The relative viscosity in this regime is a nonlinear function of bubble volume fraction with a sharp increase close to $ψ_{t}$ . The shear-induced microstructural rearrangements occurring before the threshold volume fraction affect relative viscosity to a great extent (Vermant et al., 1998, Herzhaft, 2002, Leal-Calderon et al., 2007b, Wang et al., 2015). Fig. 1 schematically shows the transition of monomodal bubble configuration that occurs in the concentrated regime.

For bubble volume fraction larger than that of the threshold packing $ψ > ψ_{t}$ , the terminology to describe the mixture changes from ”energized fluids” to ”foams”. The most important difference between energized fluids and foams happens while passing the wet limit (see Fig. 1), when the shape of bubbles changes from spherical to polyhedral (Van der Net et al., 2006, Fortuna et al., 2012) due to geometrical constraints. Thereafter, the energy at the interface between bubbles increases above the minimum surface energy. The interface energy includes the chemical energy of the substances inside bubbles and thin films (made by ambient fluid), and the energy associated with surface tension (Aubert et al., 1986, Langevin, 1999). The wet limit is defined as the state where the maximum possible volume of dispersed phase presents in fully spherical shapes that are closely packed whether hexagonally or randomly, thus $0.64 < ψ_{t} \leq 0.74$ . Note that the mentioned delineation for wet limit (or the threshold packing) is only valid for monomodal dispersion of bubbles, and as modality increases, i.e. broader bubble size distribution, this range also shifts to higher volume fraction values (Drenckhan and Langevin, 2010).

Any addition of the dispersed phase volume fraction beyond the wet limit results in osmotic pressure (i.e. tendency to take up ambient phase into the bubble configuration) (Yazhgur et al., 2016, Hhler et al., 2008) and distortion of the sphericity of bubbles to polyhedra (Aubert et al., 1986, Cohen-Addad et al., 2013). The change in the surface area of interfaces comes with an energy cost, because interfaces are in a higher energy state (i.e. molecules are in a thermodynamically unfavorable state) due to surface tension. Further increase of the volume of dispersed phase (to create foam) requires further expenditure of energy which is not tolerable for a pure continuous fluid that needs to be transformed to a network of thin films (see Fig. 1) with extremely large surface area. One way to overcome this issue and form stable foams beyond the wet limit is to add an impurity (e.g. additional surfactants) to the continuous phase. Surfactant molecules accumulate at interfaces which is an energetically favorable place for them (Aubert et al., 1986). The positioning of surfactant then reduces the surface tension (or surface energy) to some extent that improves the stability of films possessing high superficial area.

Ultimately, even in the presence of enough surfactant, increasing the volume fraction of bubbles decreases the amount of fluid available to sustain the films network, and the foam finally collapses. The onset of collapsing is referred to as dry limit beyond which the phase inversion occurs, i.e. the dispersed gas phase becomes continuous, and the liquid phase forms suspended droplets. The region beyond the dry limit is known as mist (Reidenbach et al., 1986, Herzhaft, 1999).

3. Energized fluids

3.1. Physical description

The presence of a cloud of gas bubbles in a liquid dramatically changes the transmission of stress by the bulk ambient fluid, and varies the rate at which the constituents of the mixture exchange momentum. According to the response of such fluids to external forces like shear force, their behavior may generally be categorized as Newtonian or non-Newtonian. For Newtonian fluids, the relation between shear stress and shear rate is linear, and the constant of proportionality is the shear dynamic viscosity. For non-Newtonian fluids showing shear thinning, shear thickening, and finite normal stress differences (Loewenberg and Hinch, 1996, Zinchenko, 2003), the relation between the shear stress and shear rate is nonlinear, which makes it difficult to define the shear dynamic viscosity. In these fluids the shear dynamic viscosity may increase or decrease with the applied shear rates (or shear stress) showing shear thickening or shear thinning behaviors, respectively. In almost all practical applications, energized fluids exhibit non-Newtonian behaviors, and mostly tend to shear thin. The rheological behavior of energized fluids is heterogeneous microscopically, and this is rooted in their internal microstructure (e.g. bubble size, spatial and orientation distribution) evolution under different shear conditions. For example, when flowing through a pipe where the shear stress field is inhomogeneous, the viscosity varies spatially. It has been shown that the macroscopic rheology of such fluids depends strongly on the dispersed phase volume fraction, shape, state of dispersion (Faroughi and Huber, 2017), shear history and finally the rate of deformation (Choi and Schowalter, 1975, Rust and Manga, 2002, Stein and Spera, 2002, Derkach, 2010, Faroughi and Huber, 2015, Tasaka et al., 2015). Due to the non-Newtonian nature of energized fluids, the coefficient of viscosity is measured at specific values of the shear rate to represent the so-called apparent or effective (non-Newtonian) viscosity.

To characterize the dynamics of these fluids, several dimensionless numberscan be defined that carry important information about the fluid's microstructures, (e.g. bubble size), shear conditions, thermal effect and the relative density between constituent phases. One may define the bubble Reynolds number as,(1) $R e = \frac{Inertial Forces}{Viscous Forces} = \frac{ρ_{m} R_{b}^{2} \dot{γ}}{μ_{m}},$ where $ρ_{m}$ and $μ_{m}$ are the ambient fluid density and shear viscosity, respectively, and $R_{b}$ is the bubble radius. The other important dimensionless number is the Peclet number,(2) $P e = \frac{Hydrodynamic Forces}{Brownian Force} = \frac{6 π μ_{m} R_{b}^{3} \dot{γ}}{K_{B} T},$ where $\dot{γ}$ denotes the applied shear rates, T is the absolute temperature and $K_{B}$ is the Boltzmann constant. At low Re number ( $R e ≪ 1$ ) and high Peclet number, ( $P e ≫ 10^{3}$ ), the flow dynamics of energized fluids is mainly controlled by external body forces, bubble-bubble forces, and hydrodynamic interactionsbetween bubbles, and Brownian motion is vanishing. Indeed, at low bubble Stokes number ( $S t k ≪ 1$ ) defined as,(3) $S t k = \frac{Bubble response time}{Fluid response time} = \frac{ρ_{b}}{ρ_{m}} R_{b} ς R e,$ the phase segregation to study the rheological response can be neglected. In Eq. (3), $ρ_{b}$ is bubble density and ς is a shape factor. This assumption assures that bubbles follow the bulk fluid flow streamlines.

Non-interacting gas bubbles may deform due to shearing. The spherical shape, which possesses the minimum surface energy, is deformed by a surface stress proportional to the viscous stress, $μ_{m} \dot{γ}$ . The stress that acts on the surface of bubble to restore its spherical shape is proportional to the surface tension between constituents, $σ R_{b}$ , where σ represents the surface tension. One can derive a dimensionless number between these stresses acting on the surface of bubbles to determine the equilibrium shape of bubbles, the so-called capillary number which is defined as,(4) $C a = \frac{Viscous forces}{Surface tension forces} = \frac{μ_{m} \dot{γ} R_{b}}{σ} .$

As the capillary number increases, bubbles start to deform. At early stage of deformation, i.e. at low and intermediate capillary number, the isotropy of emulsions is retained. However, at high capillary numbers, all bubbles are deformed and aligned with the imposed shear direction, which results in growing anisotropic microstructures in the fluid (Chaffey and Brenner, 1967, Manga and Loewenberg, 2001, Derkach, 2010). Deformation of a bubble to the second order in a shear flow is characterized by(5) $r - R_{b} (1 + f (\dot{γ}, r) D_{R} + f^{'} (\dot{γ}, r) D_{R}^{2}) = 0$ where r is the magnitude of the position vector from the center of a bubble to the interface between bubble at the ambient fluid, and f and $f'$ are complex functions that depend on the shear rate tensor and position vector from the center of bubbles (see Faroughi & Huber (Faroughi and Huber, 2015) and Greco (2002) for more details about these functions). In Eq. (5), $D_{R}$ characterizes the deformation of bubbles which can be obtained as (considering the viscosity ratio between constituents, $λ = μ_{d} / μ_{m}$ is very small, $λ \to 0$ ),(6) $D_{R} = \frac{R_{m a x} - R_{m i n}}{R_{m a x} + R_{m i n}} ≃ C a$ where $R_{m a x}$ and $R_{m i n}$ are respectively the maximum and minimum axis lengths of an ellipsoid obtained from a spherical bubble under deformation.

3.2. Effective viscosity of energized fluids

It is common practice to treat an energized fluid as a homogeneous Newtonian fluid under specific shearing conditions with corresponding effective viscosity. The viscosity of the homogenized or equivalent fluid is generally expected to be higher than that of the ambient fluid as the rate of energy dissipation per unit volume increases. However, if the surface tension between constituents of the energized fluids is not high enough to restore the spherical shape of bubbles under deformation, the effective viscosity may fall below than the viscosity of the ambient fluid (Pal, 2004, Stein and Spera, 2002, Manga and Loewenberg, 2001, Rust and Manga, 2002, Faroughi and Huber, 2015).

The ”homogenization” process assumes the equivalence of the work performed by the deviatoric component of the stress tensor (symbolized by $τ$ ) on the boundary of the actual energized fluids and the homogenized fluid.

Thus, the homogenization process reads as,(7) $\int (τ_{m} + N τ_{b}) u \cdot n d A = \int (τ_{ψ}) u \cdot n d A$ where $τ_{m}$ , $τ_{b}$ and $τ_{ψ}$ are respectively the ambient fluid, bubble and equivalent fluid stress tensors, $n$ is a outward unit vector normal to the surface of bubble and $u$ is the fluid flow at large distances from the bubble center. The deviatoric stress tensor for the energized fluid (left-hand-side of Eq. (7)) consists of the stress associated with ambient fluid and the stress associated with the disturbance due to the presence of a bubble summed over the number of bubbles, N. The latter part of the stress tensor may include different orders of bubble deformation. The zeroth order of deformation ( $C a ≪ 1$ ) was calculated by Taylor (1932), and the first and the second order of deformation was estimated in later studies by several authors (Chaffey and Brenner, 1967, Schowalter et al., 1968, Greco, 2002).

3.2.1. Dilute energized fluids

The main assumption embedded in Eq. (7) (i.e. the linear summation of the effect of bubbles on the bulk stress field) is strictly valid for dilute regimes where bubbles do not interact at all. The macroscopic relative viscosity of energized fluids predicted using this assumption is linearly proportional to the bubble volume fraction, especially when $ψ \to 0$ . For semi-dilute regimes, the non-linear effect of bubble perturbation of the bulk stress field was considered by Batchelor & Green (Batchelor and Green, 1972), who found that the relative viscosity becomes a quadratic function of the bubble volume fraction. At higher volume fractions of bubbles, the strong hydrodynamic and inter-bubble (e.g. lubrication, collision and friction) interactions need to be considered within the stress disturbance function. Theoretical approaches are unable to describe these highly complex many-body interactions, and thus, experiments (Pal, 2001, Pal, 2003, Rust and Manga, 2002, Stein and Spera, 2002, Murai et al., 2015) and numerical studies (Kang et al., 1997, Zinchenko, 2003, Tran-Duc et al., 2013, Wang et al., 2015) are utilized to constrain the highly non-linear relative viscosity of intermediate and concentrated energized fluids. Phenomenological approaches have also been used widely to characterize concentrated systems. See Table 1 for a non-exhaustive list of well-known relative viscosity models for concentrated energized fluids and their range of applicability based on the volume fraction of dispersed phase as well as deformation parameters.

Table 1. A non-exhaustive list of well-known viscosity models to predict the relative viscosity of energized fluids. The range of bubble volume fraction, ψ, and capillary number,

C a

, over which they intend to be applied is also reported (Faroughi and Huber, 2015).

Reference	Equation	Specific Parameter	Concentration	Deformation
Taylor (Taylor, 1932)	$\frac{μ_{ψ}}{μ_{m}} = 1 + ψ$	–	$ψ \to 0$	$C a ≪ 1$
Mackenzie (Mackenzie, 1950)	$\frac{μ_{ψ}}{μ_{m}} = 1 - \frac{5}{3} ψ$	–	$ψ \to 0$	$C a > 1$
Oldroyd (Oldroyd, 1959)	$\frac{μ_{ψ}}{μ_{m}} = \frac{1 + \frac{3}{5} ψ}{1 - \frac{2}{5} ψ} [\frac{1 + K_{1} K_{2} {(\frac{6 C a}{5})}^{2}}{1 + K_{1}^{2} {(\frac{6 C a}{5})}^{2}}]$	${\begin{matrix} K_{1} = \frac{1 + \frac{2}{3} ψ}{1 - \frac{2}{5} ψ} \\ K_{2} = \frac{1 - ψ}{1 + \frac{3}{5} ψ} \end{matrix}$	Up to 20%	All range
Frankel & Acrivos (Frankel and Acrivos, 1970)	$\frac{μ_{ψ}}{μ_{m}} = 1 + ψ [\frac{1 - \frac{12}{5} C a^{2}}{1 + {(\frac{6}{5} C a)}^{2}}]$	–	Up to 10%	All range
Ducamp & Raj (Ducamp and Raj, 1989)	$\frac{μ_{ψ}}{μ_{m}} = e x p (- k \frac{ψ}{1 - ψ})$	$2.5 \leq k \leq 4$	Up to 45%	$C a > 1$
Pal (Pal, 2004)	$\frac{μ_{ψ}}{μ_{m}} = {(1 - \frac{ψ}{ψ_{t}})}^{- ψ_{t}} (\frac{1 + K_{1} K_{2} {(\frac{6 C a}{5})}^{2}}{1 + K_{1}^{2} {(\frac{6 C a}{5})}^{2}})$	${\begin{matrix} K_{1} = {(1 - \frac{ψ}{ψ_{t}})}^{\frac{- 16}{15} ψ_{t}} \\ K_{2} = {(1 - \frac{ψ}{ψ_{t}})}^{\frac{8}{5} ψ_{t}} \end{matrix}$	Up to $ψ_{t}$	All range
Llewellin & Manga (Llewellin and Manga, 2005)	$\frac{μ_{ψ}}{μ_{m}} = {\begin{matrix} {(1 - ψ)}^{- 1} \\ {(1 - ψ)}^{5 / 3} \end{matrix}$	Minimum model	Up to 7%	$\begin{matrix} C a \leq 1 \\ C a > 1 \end{matrix}$
	$\frac{μ_{ψ}}{μ_{m}} = {\begin{matrix} 1 + 9 ψ \\ {(1 + 22.4 ψ)}^{- 1} \end{matrix}$	Maximum model	$\begin{matrix} Up to 50 % \\ Up to ψ_{t} \end{matrix}$	$\begin{matrix} C a \leq 1 \\ C a > 1 \end{matrix}$
Faroughi & Huber (Faroughi and Huber, 2015)	$\begin{matrix} μ_{r} {(\frac{1 - \frac{12}{5} C a^{2} {(μ_{r})}^{2}}{1 - \frac{12}{5} C a^{2}})}^{- \frac{4}{5}} = \\ {(1 - Ω \frac{ψ}{1 - ψ})}^{- Ω^{- 1}} \end{matrix}$	$Ω = \frac{1 - ψ_{t}}{ψ_{t}}$	Up to $ψ_{t}$	All range