Advances in the Rheological Characterization of Slurries of Elongated Particles †

Inherent challenges regarding the rheological characterization of slurries of elongated particles have necessitated the development of alternatives to standardized rheometers. These methods of measurement, and the associated advances in the quantification of shear and normal stress measurements, are described. Also, recent advances in modeling and predictive capabilities are summarized. During shearing flows, confinement substantially influences the orientation distribution of the particles; this change in the microstructure impacts the rheology, even as the smallest confining dimension exceeds seven particle lengths. The slow development of the orientation distributions renders additional difficulties in evaluating the rheology. Achievements of the measurement methods include a universal shear viscosity as a function of concentration for a wide range of particle lengths to diameters (aspect ratios). The jamming limit (divergence of the viscosity with concentration) of the suspensions has been also shown to scale differently than for spheres. More general dynamics of the suspensions and the additional needs for measurement improvements are discussed.


Introduction
Suspensions, a class of complex fluids wherein an insoluble phase, either a solid or a liquid is dispersed in a liquid (Larson, 1999), are ubiquitous in nature. Suspended particles can be spherical or near-spherical, discs, or elongated in the form of fibers, and can be rigid or soft and flexible. Suspensions occur naturally, with macroscale examples such as lava and silt-laden rivers where the solid particles are anisotropic and polydisperse. At the microscale, blood is an example of another natural suspension where disc-like platelets are suspended in plasma. These suspensions can be further categorized on the basis of their particle sizes as Brownian and colloidal with sizes on the order of tens of nanometers to a few microns, to non-Brownian and non-colloidal for particle sizes larger than a few microns.
In industry, suspensions of non-colloidal particles are found in a wide variety of applications. Transport and processing of slurries containing wood fibers in the manufacture of paper products is an example where the particles are long and flexible (Lundell et al., 2011). Rod-like particles and flexible fibers are added to materials to create compos-ite materials with enhanced strength and performance, as with fiber-reinforced concrete (Hassanpoura et al., 2012). Fibers, both rigid and flexible, are also added to modify the rheological properties of fluids used in the production of oil (Bivens et al., 2005;Elgaddafi et al., 2012;Hammond, 1995;Osiptsov, 2017). In contrast, the presence of fibrous minerals can have a detrimental impact on the operation of flotation processes (Patra et al., 2012;Somasundaran et al., 2019), which are often used at mining sites to separate high value minerals from a slurry.
Various trade-offs exist while designing and operating processes that contain particle suspensions. Decisions often center on the balance between high energy usage at high particle concentrations and high suspending fluid usage at low particle concentrations. Regarding the latter, concerns regarding overuse of water resources have been associated with many industries, including the example of mining and the operation of on-site flotation processes (Liu et al., 2011;Rao and Finch, 1989). At the other limit of high solids content, the viscosity of the suspension increases and even diverges as the solids content approaches maximum packing. As an immediate consequence, the power requirements for pumping and mixing suspensions can increase significantly, thereby significantly increasing processing costs.
Reliable, quantitative information on the rheology of suspensions can assist in managing these trade-offs. For suspensions of spheres that are non-Brownian and noncolloidal, there has been considerable development of rheological models. For volume fractions ϕ < 0.02 (i.e the the ratio of the total volume of the suspended particles to the total volume of the suspension), the viscosity η s increases linearly according to Einstein's equation (Einstein, 1906), where η f is the viscosity of the suspending fluid. Higher order corrections for suspension viscosity have been determined which extend the range of the volume fraction (Batchelor and Green, 1972). However, the complexity of the problem quickly surpasses analytical capabilities as the concentration is increased, and computational methods such as the Stokesian dynamics method (Brady and Bossis, 1988) must be employed in order to make progress on calculations of rheology based on first principles. Experimentally characterizing the rheology of concentrated suspensions of spheres can pose difficulties as well. Acrivos (1995) summarized discrepancies among published reports on the shear viscosity and attributed the differences to a number of effects including shear-induced migration (Leighton and Acrivos, 1987), which causes a non-uniform distribution of particles within the testing geometry. One generally accepted correlation for the shear viscosity of suspensions of non-Brownian spheres in Newtonian fluids over a wide range of volume fractions is given by the Krieger-Dougherty equation (Krieger and Dougherty, 1959), where ϕ m is the maximum volume fraction at which the particles pack randomly and [η] ≈ 2.5. Still, some debate persists regarding rheological properties of concentrated suspensions of spheres, such as the sign of the first normal stress difference (Guazzelli and Pouliquen, 2018). Making well-informed decisions regarding processes that utilize suspensions of particles that are non-spherical is hampered by a relative lack of data and incomplete information on how to interpret that data. In the specific case of suspensions of elongated particles, though the viscosity scales with an increasing particle volume fraction, the relationship is much more complex. As an example, Fig. 1 shows rheological data obtained experimentally by Tapia et al. (2017) for rigid polyamide fibers of different aspect ratios and particle diameters. It is evident that the particle shape has a significant impact on the rheology, and hence equipment design as well as operating decisions must consider these factors. However, obtaining reliable and accurate rheological data can be challenging for these types of suspensions composed of large rigid fibers at high concentration owing to confinement effects and the tendency of the suspension microstructure to change drastically depending upon the pattern of shearing flow.
This review focuses on the concentrated regime of slurries of elongated particles in effort to consolidate the sparse work done previously in this area and to summarize recent and ongoing studies. Section 2 outlines the properties of the particles and fluids, as well as the regimes of particle concentrations in the slurries under consideration in this review. Section 3 examines the rheological techniques and measurements in the semi-dilute (Section 1) and concentrated (Section 2) regimes of concentration, focusing on studies of the near-jamming limit. This section also explores the effects of time-varying flow fields, focusing on the relatively recent advances in oscillatory flows (Section 2). Normal stress differences, which are of importance in predicting shear-induced migration, are discussed in Section 4. Some comments on remaining challenges are offered in the conclusions.

Suspension properties
This review focuses on slurries of elongated particles in viscous Newtonian fluids. Many of the studies described use rod-shaped particles of aspect ratio, A = L/d >> 1, where L is the length of the rod and d is its diameter. Unless specified otherwise, these particles are sufficiently large to be non-Brownian. Consequently, the rotational Péclet number Pe r = γ is the rate of shear and d r is the rotational diffusivity, for a rod of high aspect ratio. The fluid viscosity is η f , k is the Boltzmann constant, and T is the temperature. An  Tapia et al. (2017) showing the dependence of the shear viscosity η s on the particle volume fraction ϕ for rigid fibers. The measured viscosities depend on the particle aspect ratio A = L/d, where L and d are the particle length and diameter, respectively. The results indicate that the viscosity diverges at a different volume fraction for each aspect ratio, which roughly corresponds to maximum random packing (ϕ m ).
alternative definition for the rods to be non-Brownian requires a large Péclet number as defined by Pe = L 2 γ (31) /D T , where D T is the mean center of mass diffusivity. Basing the definition of non-Brownian on Pe >> 1 is more stringent than using Pe r since Pe = 9Pe r .
The particles chosen for idealized experiments are monodisperse in length and diameter. The materials of these fibers are generally polymeric, for ease of production, use, and manipulation. Polyamides like Nylon-6 and polystyrene (PS) are used in conventional rheometry experiments owing to the ease of obtaining particles with a wide array of lengths, diameters, and aspect ratios, and mechanical strength. For imaging experiments, poly-methylmethacrylate (PMMA/acrylic) obtained from the core of fiber optic cables have been optimal owing to their uniform optical properties.
The suspending media are Newtonian in all of the work described in this review. The density of the fluids (or mixtures) were selected or altered to match the density of the particles; the particles neither sediment to the bottom or float to the top, negating any possible effects of gravity on the dynamics and rheology. The above-mentioned conditions are imposed upon the slurries to simplify the subsequent measurements and analyses, by assuming that in the absence of Brownian fluctuations and gravity, the motion of the particles in the slurries are primarily driven by the flow imposed upon them. For imaging experiments, the suspending media must also be chosen, or its composition altered, in order to match the refractive index of the particles. When the fluid is subsequently dyed with a fluorescing agent and properly illuminated, the contrast can be used to obtain spatial and orientation data.
The elongated particles on which this review focuses are rigid, neither stretching nor bending under the forces imposed by the flow. For materials that are relatively soft, shearing forces may be sufficient to bend a fiber and impact the dynamics and rheology. Consequently, criteria have been developed for estimating the stress at which fibers bend appreciably. Forgacs and Mason (1959) determined a critical stress, Σ crit , at which a slender rod buckles under compressional flow as, where E b is the bending modulus, estimated to be twice the Young's modulus of the particle. This model provided good agreement for the bending of polymer fibers in shearing flows with stresses of η f γ (31) (Forgacs and Mason, 1959). It is assumed that any deformation of the elongated fiber is negligible when the maximum stress exerted by the flow field is much smaller than Σ crit . A number of authors have developed similar ideas, defining a ratio of the imposed fluid stresses of η f γ (31) and critical stresses as calculated using relations identical or similar to Eq. (4). Switzer and Klingenberg (1959) defined an effective stiffness, S eff , as where E y is the Young's modulus of the particle. Others (Lauga, 2007;Wandersman et al., 2010) inverted the ratio of shear stress and critical stress to create the "sperm number" S P as, Young and Shelley (2007) Despite the differences in the various definitions, bending of a fiber can be ignored so long as the product of the imposed shear and viscosity (η f γ ) are much smaller than the stress at which the fiber buckles.
The rheology and dynamics of suspensions of rigid fibers differ greatly depending on the extent to which rotation is impeded by the presence of other fibers. Consequently, (Doi and Edwards, 1978) suggested defining concentration regimes in terms of simple geometric constraints; examples are shown in Fig. 2 of fibers randomly placed in a fixed volume at dilute, semi-dilute, and concentrated conditions.

Fig. 2
Images of rigid fibers of aspect ratio A = L/d = 10 placed randomly in a volume at concentrations within the a) dilute (nL 3 << 1), b) semi-dilute (1/L 3 < n < 1/L 2 d), and c) concentrated (nL 2 d > 1) regime. The volume fraction ϕ corresponding to each dimensionless number density is given in the figure. At the highest concentration of nL 2 d = 3.0, the suspension is isotropic for this aspect ratio.
For the case of dilute suspensions (nL 3 << 1), spherical volumes that enclose each particle do not generally overlap and each fiber can rotate freely. Here n = N/V is the number density, the number of particles N per unit volume V of the suspension. Note the dimensionless number density nL 3 is equivalent to 4A 2 ϕ/π, where ϕ is the particle volume fraction. In the semi-dilute region, fibers can no longer rotate freely. The upper limit of the semi-dilute regime (nL 2 d = 4Aϕ/π ≈ 1) roughly corresponds to close packing of a collection of flattened discs having a diameter of L and thickness d. Further increases in concentration beyond nL 2 d = 1 can occur while the spatial and orientation distribution remains isotropic, as demonstrated in Fig. 2c for A = 10 and nL 2 d = 3.0. However at some point, maximum random packing is exceeded and adding more particles requires organizing the particle orientations into a nematic phase. Studies (Williams and Philipse, 2003) indicate that the maximum volume fraction at which sphero-cylinders can still pack randomly scales inversely with the aspect ratio. For the case of A = 10 shown in Fig. 2, maximum random packing occurs at nL 2 d ≈ 4.7, or ϕ ≈ 0.37; higher concentrations requires aligning the fibers. At A = 20, the transition occurs at ϕ ≈ 0.22, or nL 2 d = 5.6.

Shear rheology and structure
The rheology of suspensions of well-dispersed rigid fibers depends strongly upon the orientation distribution. As a quantitative example, the shear viscosity for dilute solutions of rigid fibers (nL 3 < 1) that are force-free (i.e. non-Brownian) is given by Batchelor (1970)   3 2 2 s f 1 , where η s is the viscosity of the suspension relative to that of the fluid (η f ), c is a constant that depends on the particle size and shape, and the expression ignores the small contributions to stress from the thickness of the particle, which is assumed to be of high aspect ratio. The brackets,  , indicate an average over all fibers in the suspension, where the fibers have an orientation of p. Here, p x and p y represent the orientation components in the flow (x) and gradient (y) directions.
Generally, the orientation (and position) of each rod in a suspension is set by balancing the non-hydrodynamic torques (and forces) with hydrodynamic drag. In the particular case of a force-free fiber at dilute conditions, the fiber tumbles indefinitely in a well-defined orbit that depends on its initial orientation (Jeffery, 1922). Consequently the viscosity η s oscillates in time as indicated in Eq. (8) owing to its dependence of the ever-changing orientation moment . Measurements performed by Ivanov et al. (1982) confirm that the viscosity oscillates during steady shear, though the amplitude of the oscillation decays over time. The decay to a steady value is likely due to small effects that influence the torque balance on the fibers, possibly including weak Brownian motion (Hinch and Leal, 1972) or inertia (Einarsson et al., 2015).
Thorough reviews of the rheology of dilute suspensions are available (Powell, 1991;Petrie, 1999) that discuss these issues in depth. The remainder of this section describes the measurement of more concentrated suspension where particle interactions begin to influence the orientation distributions and the measured rheology. Fig. 3 shows viscosity measurements by Bibbó (1987) and Djalili-Moghaddam and Toll (2006) for glass and nylon fibers of aspect ratios A ≈ 32 and 50. Measurements are shown for volume fractions of ϕ = 0.01 to 0.1, which correspond to dimensionless number densities nL 3 spanning the semi-dilute concentration regime (nL 3 > 1 and nL 2 d < 1).

Fig. 3
Shear viscosity η s relative to the fluid viscosity η f as a function of fiber concentration nL 3 for aspect ratios of a) A ≈ 50 and b) A ≈ 32; the concentration range spans the semi-dilute regime of nL 3 > 1 and nL 2 d < 1. Experimental measurements of Bibbó (1987) and Bounoua et al. (2016a) are compared to simulations (Mackaplow and Shaqfeh, 1996;Claeys and Brady, 1993) and the theory for semi-dilute suspensions of rigid rods by Shaqfeh and Fredrickson (1990). A value of 2 2 x y p p [51] = 0.01 was used when evaluating the theory numerically. Error bars are shown for the sets of data when available; for the data of Mackaplow and Shaqfeh (1996), the error is often smaller than the size of the data point.
Viscosities for higher aspect ratio fibers are expected to be larger at the same concentration and for the same orientation distribution, as indicated in the figure by the semi-dilute theory (Shaqfeh and Fredrickson, 1990) where the same value of 0.01 for 2 2 x y p p (32) was used in calculating the curves. Indeed, the experimentally measured viscosities in Fig. 3 are larger for the higher aspect ratio fibers. However, the influence of the orientation distribution is strong: the viscosity of suspensions with fibers aligned with the flow was significantly lower than that of an isotropically distributed fibers as seen in Fig. 3a (isotropic versus non-isotropic data by Bibbó (1987)). For this reason, the influence of bounding walls in confined geometries must be considered when evaluating the measured rheology, especially as the concentration approaches and surpasses nL 2 d = 1.
Beyond the dilute limit of nL 3 = 1, interactions between particles become increasingly important, and predicting the rheology depends on at least two complicated calculations. First, accurate determination of the spatial and orientation distribution is necessary. Second, calculations that relate the microstructure and stresses must be performed. In both cases as reviewed in detail by Butler and Snook (2018), hydrodynamic interactions between the particles must be considered. This differs from the dilute and concentrated regimes where the particles are either too far apart or the dynamics and rheology are dominated by direct particleparticle contacts. Calculations that consider hydrodynamic interactions between rigid rods, such as those shown in Fig. 3, are computationally expensive, hence predicting the microstructure and rheology for large numbers of fibers in the semi-dilute regime remains an outstanding challenge.
Some results shown in Fig. 3 represent the high shear limit of the measured viscosity. In the semi-dilute regime, many slurries of rigid fibers exhibit shear-thinning rheology, where the measured viscosity decreases with increases in rate of shear (Ganani and Powell, 1985;Powell, 1991;Bricker et al., 2008;Bounoua et al., 2019). One example from Bricker et al. (2008) is shown in Fig. 4 for a suspension of ellipsoids with a moderate aspect ratio of A = 7 at a concentration of nL 3 = 4.4 (ϕ = 0.07). Although the particles were relatively small in these experiments, the value of the rotational Péclet number was large, ranging from approximately Pe r = 10 3 at the lowest shear rate of γ (31) = 0.1 s -1 to 10 6 at 100 s -1 . Note that Fig. 4 also shows viscosity data over the same range of shear rates for spheres having properties (material of composition, volume of each particle, and method of manufacture) identical to that of the ellipsoids.
The substantial shear-thinning of the suspension of ellipsoids versus the spheres indicates a remarkable difference in the rheology originating from the shape of the particles. These results contradict numerical simulations (Dinh and Armstrong, 1984;Shaqfeh and Fredrickson, 1990;Schiek and Shaqfeh, 1995) which predict viscosities that are independent of shear-rate. Indeed, the rheology of suspensions of non-Brownian and non-colloidal particles is expected to be pseudo-Newtonian regardless of their shape: though attaining steady-state may require a large strain, or time, the steady viscosity of the fluids should be independent of shear rate. Many ideas regarding the origin of the shear thinning have been given, including the possibility that non-Newtonian suspending fluids are responsible (Bibbó et al., 1985;Maschmeyer and Hill, 1977) and that the fibers used in the experiments deform under the imposed stresses (Sepehr et al., 2004;Bennington et al., 1990). Attractive colloidal forces have also been incorporated into models to explain shear thinning (Chaouche and Koch, 2001;Bounoua et al., 2016b), and this line of investigation has been argued to be promising (Butler and Snook, 2018).

Concentrated suspensions
The rheology of suspensions of rod-shaped elongated particles in a Newtonian fluid beyond the semi-dilute regime remain largely unexplored. This is largely due to the difficulty in performing rheological experiments with large, rigid fibers at dimensionless number densities nL 2 d above 1.0. Under these conditions, confinement is a significant issue. Standard rheological tooling employs shearing geometries having gaps between the bounding walls of H≲1 mm. This is particularly problematic for suspensions composed of non-colloidal fibers where the particles of interest can have lengths L ≳ 0.1 mm. Typically it has been assumed that the rheological measurements are not influenced by the bounding walls so long as H/L > 2 or 3 (Bibbó et al., 1985;Powell, 1991;Chaouche and Koch, 2001). This is true for dilute, or even semi-dilute, suspensions since the additional drag on any one rod rapidly decays as the distance from the bounding wall increases (Zurita-Gotor et al., 2007;Park and Butler, 2009). However, at higher concentrations the Fig. 4 Steady shear viscosity η s as a function of shear rate for a suspension of ellipsoids of aspect ratio A = 7 and for spheres (Bricker et al., 2008). Both suspensions have a volume fraction of ϕ = 0.07, and the particles differ only in their shape, as each particle was composed of an identical volume of polystyrene and was manufactured using the same process.
bounding walls have a strong effect on the rheology: the walls prevent particles from rotating and moving freely, altering the orientation distribution. This disturbance in the structure imposed by the wall propagates a significant distance into the fluid at high concentrations. The measured rheology will not match that of a bulk fluid unless the wall spacing is much larger than 3L.
Consequently, Tapia et al. (2017) used a large, custombuilt device to make the viscosity measurements shown in Fig. 1. The rheometer was originally constructed by Boyer et al. (2011a) and was used by them and Dagois-Bohy et al. (2015) to characterize a Newtonian fluid containing spheres. This "imposed load" rheometer consisted of a wide-gap annular shear cell with a movable top plate, as shown in Fig. 5. The movable top plate has a porous screen which allowed fluid to pass freely across the plate, but prevented fibers from passing. This enabled the dynamic adjustment of the particle volume fraction while maintaining the volume of the particles.
Experiments were conducted in two modes: conventional rheometry with a controlled volume and a novel "pressure-imposed" mode. For the former mode, the volume fraction of the suspension was maintained at a constant value by holding the top plate at a fixed height; the torque was measured as a function of rotation rate while also recording the normal force needed to prevent the top plate from moving. In the pressure-imposed mode, a constant normal force was applied to the top plate, allowing it to move while measuring the torque as a function of rotation rate. In this case, the volume fraction ϕ varied with the rate of rotation. The torque and rotation rate were used to calculate the shear viscosity η s as a function of ϕ as shown in Fig. 1. The normal force on the plate was corrected for buoyancy and divided by the plate area to give a pressure P. This quantity represents the particle pressure that the fluctuating particles exert on the porous plate and is reported as a normal viscosity, as was done by Morris and Boulay (1999). These results are shown in Fig. 6. A unique result provided by these experiments is the measurement of the rheology at concentration ranges approaching maximum packing fraction, ϕ m . These measurements were made possible by performing experiments in the "pressure-imposed" mode. Fig. 1 and Fig. 6 show that the shear and normal viscosities diverge at a value of ϕ that depends on the aspect ratio A of the fibers and corresponds closely to maximum random packing (ϕ m ) for each aspect ratio. As the concentration, ϕ approaches ϕ m , the suspension approaches the "jamming transition". The maximum   (Tapia et al., 2017). The normal viscosity is a unique measurement, as it is a bulk rheological measure of particle pressure.
packing fraction serves as an important scaling factor as well: as shown in Fig. 7, the rheological data collapses onto single-valued curves over the entirety of the range of particle sizes, aspect ratios, and concentrations measured in the experiments. Constitutive equations (Tapia et al., 2017) for the shear and normal viscosities, along with the coefficient of friction, developed using the maximum packing fraction as a scaling factor, are where μ s is the coefficient of static friction, c 1 and c 2 are constants obtained from curve-fitting, and differ from the values obtained from experiments with spherical particles (Dagois-Bohy et al., 2015). Shear-thinning rheology has also been reported for suspensions of elongated fibers in the concentrated regime and there is a wider range of measured values in this range of concentrations (Ganani and Powell, 1985;Powell, 1991). Theories similar to those mentioned at the end of Sec. 3.1 have been put forward to explain both observations, including the presence of attractive interparticle forces. Tapia et al. (2017) did not report shear-thinning behavior in their system, but did report finite yield stresses as have others. These yield stresses increased with increasing volume fraction and fiber aspect ratio. Many have attributed the yield stresses to attractive forces between the fibers, just as with shear thinning (Mongruel and Cloitre, 1999;Chaouche and Koch, 2001;Bounoua et al., 2016b).

Unsteady shear rheology
As with other non-Newtonian fluids such as polymer and colloid solutions, rheological characterizations of fiber suspensions under time-varying shearing flows can probe the microstructure and forces that act within the system. However, the absence of Brownian motion and a welldefined equilibrium state for these fluids requires careful consideration regarding what is measured; for example, the Cox-Merz rule can not be relied-upon for interpreting small amplitude oscillatory flows. Here, the shear reversal and oscillatory shear rheology for concentrated suspensions of fibers are described. As with steady shearing flows, a point of concern is the effect of confinement on the measurements. Bounoua et al. (2019) measured the time-dependent shear viscosities upon suddenly reversing the direction of a steady shearing flow. Fig. 8 shows the measurements as a function of strain (γ = γ (31) t) for a concentrated suspension (nL 2 d = 3.2, ϕ = 0.14) of rigid fibers of aspect ratio A = 18;   10) and (11). The insets are log-log plots to determine the exponent value in Eq. (10) and is compared to that of spherical shaped particles (dashed).

Shear reversal
the results share similarities with other studies of shear reversal of fiber suspensions (Sepehr et al., 2004). First, the samples were sheared until the viscosity reached steady state. Note that these suspensions of fibers were found to be shear thinning, consequently the steady shear viscosity differed depending upon the applied stress. The reversal was initiated at γ = 0 (i.e. t = 0), at which point the measured viscosity dropped significantly. The low point following the reversal was followed by a transient overshoot in the viscosity, before it returned to its steady value as confirmed by comparing viscosity values at large γ to the values at γ < 0, before reversing the flow direction.
Suspensions of non-colloidal spheres also exhibit reduction in the shear viscosity when the flow direction is suddenly changed following a period of steady shear (Gadala-Maria and Acrivos, 1980;Narumi et al., 2002;Blanc et al., 2011). Stresses resulting from the shearing flow of non-Brownian and non-colloidal particle suspensions have two sources: hydrodynamic and contact forces. The hydrodynamic contribution from each particle arises from the difference in velocity of the particle with the imposed shear flow and velocity disturbances caused by the other particles, which depends upon the relative arrangement of the particles. Direct contact forces between particles occur due to the breakdown of lubrication interactions, possibly because of surface roughness (Rampall et al., 1997;Ingber et al., 2006;Pham et al., 2015).
For shear-reversal of suspensions of spherical particles, detailed analysis (Peters et al., 2016) of Stokesian dynamics simulations indicate that the viscosity measured immediately after the reversal represents the purely hydrodynamic contributions to the shear viscosity, as the contact contribution disappears nearly entirely. During steady shearing flow, contact interactions build-up among the particles and contribute to the shear stresses. Suddenly reversing the direction of flow relieves the contact forces, which are purely repulsive, and the shear stress accordingly

drops.
Similarly, the sudden drop in the viscosity shown in Fig. 8 indicates the presence of contact interactions between the fibers. However, the overshoot in viscosity represents an important difference in the results as compared with spheres, where the viscosity returns to its steady value in a monotonic fashion. Bounoua et al. (2019) attributed the overshoot to weak attractive force between the particles, rather than the dynamics changes in the orientation distribution of the particles which surprisingly changes little during the course of the reversal process (Cieslinski et al., 2016). Measurements on suspensions of fibers that do not shear-thin (i.e. have no attractive forces) would presumably show no overshoot in the viscosity after shear reversal.

Oscillatory shearing
A suspension of rigid fibers undergoing oscillatory shear may yield vastly different measurements of the rheology, depending on the concentration, confinement, and strain amplitude of oscillation. Fig. 9a shows a strong qualitative difference for the latter case, where η'' (i.e. imaginary component of the complex viscosity) vanishes with increasing oscillation number (N) or remains finite, depending on the strain amplitude γ 0 of the oscillation (Franceschini et al., 2014). The rheological data reported is for aspect ratio A = 11 and a concentration of ϕ = 0.20; importantly, the suspension was highly confined (H/L = 1.5).
Measurements (Franceschini et al., 2011;Strednak et al., 2021) and simulations (Snook et al., 2012;Strednak et al., 2021) of the orientation distribution under oscillatory flows also demonstrate a remarkable transition in the alignment as a function of strain amplitude. As shown in Fig. 9b for the same aspect ratio A = 11 and a concentration of ϕ = 0.20, the order parameter S θ is maximized when γ 0 ≈ 2.5 (γ 0 ϕ/Φ c ≈ 1.25 where Φ c = 0.4). Here,  (Franceschini et al., 2014). The result for γ 0 < 2.60 marks a transition between vanishing and persistent values for η''. b) Simulated predictions (Strednak et al., 2021) and an experimental measurement (Franceschini et al., 2014) of the microstructure S θ at steady state as a function of the rescaled strain amplitude, γ 0 ϕ/Φ c . A value of 0.4 was used for Φ c .
where θ is the angle between the fiber projection in the flow-vorticity plane and the flow direction as depicted in Fig. 10a; the brackets (  ) indicate an average over all fibers. Note that fibers are perfectly aligned in the flow and vorticity directions when S θ = -1 and 1, respectively; the gradient direction is perpendicular to the flow-vorticity plane and is the same direction of the confinement in the experiments. Consequently, the large value of S θ at γ 0 ≈ 2.5 indicates strong alignment of the fibers in the vorticity direction, as can be directly seen in Fig. 10b which compares the initial and final distribution of the fibers predicted by simulations. Note that the simulated values of S θ quantitatively matched the experimental results over a wide range of conditions (Strednak et al., 2021). Franceschini et al. (2014) suggested that the transition in S θ and η'' is controlled by a critical effective volume fraction (Φ c ) of the suspension. The effective volume fraction of a fiber suspension is the ratio of the volume swept out by all fibers during a cycle to the total volume of the suspension, and is therefore proportional to the strain amplitude γ 0 . The critical effective volume fraction, Φ c defines the combination of volume fraction ϕ and strain displacement γ 0 at which the particles are forced to contact. For oscillatory shear below Φ c , Franceschini et al. (2014) speculated that η'' disappears as the suspension rearranges to avoid contact, though simulations indicate that at least some contacts persist (Snook et al., 2012;Strednak et al., 2021). Above Φ c , contacts between particles must occur during the course of every oscillation and η'' remains finite. As seen in Fig. 9a for the case of A = 11 and ϕ = 0.20, the transition in η'' occurs between γ 0 = 2.0 and 2.6. Consequently, Franceschini et al. (2014) defined Φ c = γ 0 ϕ ≈ 0.4. This corresponds to the value of γ 0 ϕ at which S θ is maximized for the same conditions of A = 11 and ϕ = 0.20 (blue circles in Fig. 9b). Strednak et al. (2021) found that the maximum value of S θ declined when the gap spacing H was increased while holding the fiber length L constant, both in experiments and simulations. This indicates that the vorticity alignment of the fibers is promoted by the bounding walls, as was speculated by Snook et al. (2012). This change in orientation distribution with H/L likely impacts the complex viscosity in a qualitative fashion as well. However, rheological data for oscillatory flow currently is lacking for large values of H/L.

Normal stresses
Suspensions can exhibit finite normal stress differences during shear. This non-Newtonian property of suspensions holds true even if the suspending fluid exhibits no normal stresses (i.e., the suspending fluid is Newtonian). Hence the normal stress differences originate with the particles, which must be present at a sufficiently high concentration for the normal stresses to be measurable. Non-zero normal stress differences can cause a number of effects, including shear-induced migration of particles within non-linear flows (Nott and Brady, 1994;Nott et al., 2011;Snook et al., 2016;Strednak et al., 2018).

Measuring normal stresses
Measuring normal stresses for suspensions of rigid fibers requires addressing at least two problems. As with other rheological characterizations, confinement of the particles at high concentrations creates difficulties when using standard shear cells. Additionally, the normal stresses are small and scale linearly with the rate of shear rather than quadratically, as with polymer solutions. As a result, normal stress differences can not be measured at low rates of shear and concentration using most rotational rheometers, and inertial effects can alter the measurements if the shear rate is increased too much in an attempt to generate a larger force that can be resolved by the transducers .   Fig. 10 a) The angle between the fiber's projection in the flow-vorticity plane and the flow direction is θ. The x, y, and z directions are the flow, gradient, and vorticity directions, respectively. b) Reproduced rendering of simulations from Strednak et al. (2021) showing an example initial (N = 0) and steady (N = 4000) distribution of the fibers for H = 1.5L, A = 11, ϕ = 0.20, and γ 0 = 2.5. Larger values of S θ reflect an increased alignment of fibers in the vorticity direction.
Two approaches for overcoming the challenges of experimentally measuring the normal stresses have been employed. The first makes use of the fact that normal stress differences deform the interface during slow free-surface flows as illustrated in Fig. 11. Because the geometries are large, confinement does not influence the measurement. Analysis and experimental approaches using free surface flows originally were developed for viscoelastic fluids (Beavers and Joseph, 1975;Wineman and Pipkin, 1966;Tanner, 1970). Also, the rotating rod geometry (Zarraga et al., 2000;Boyer et al., 2011b) and tilted trough (Couturier et al., 2011;Dai et al., 2013) previously were used to quantify normal stresses in suspensions of spheres. Whether or not the particles are spherical or rod-like, the analysis relating the interface deformations and normal stresses is similar, though also simpler than the analysis for elastic fluids. Each normal stress component is assumed to be proportional to the magnitude of the shear stress, |τ|, and the results are reported as coefficients of the normal stress differences, 1 2 1 2 and , where N 1 and N 2 are the first and second normal stress difference. Using either a domain perturbation method or a "cut-away" analysis, while accounting for the tendency of surface tension to flatten the surface, small interface deformation due to normal stress differences can be predicted accurately. The deformation of interface in the tilted trough is related to α 2 , while analysis of the interface in the rotating rod geometry yields α 2 + α 1 /2. Comparing the predicted deformations to those measured in experiments (see Fig. 11) enables determination of both coefficients.
The values of α 2 as determined from the tilted trough are shown in Fig. 12a (Snook et al., 2014). The values are strongly negative for small aspect ratios and high concentrations. In the case of the rotating rod geometry, Snook et al. (2014) reported that no deformation of the surface was observed across all tested conditions. Using the fact that α 2 + α 1 /2 ≈ 0 together with the results for α 2 gives α 1 > 0 for small aspect ratios and high concentrations, as shown in Fig. 12b. The result contrasts with suspensions of spheres, where the interface dips near the rotating rod (Zarraga et al., 2000) and α 1 is zero or weakly negative.
The second approach for measuring the bulk normal stresses of rigid fibers used a custom-built torsional flow cell (Bounoua et al., 2016a). As originally described by Dbouk et al. (2013), the device uses an unusually large sep- The value of α 2 + α 1 /2 can be determined from the deformation of the interface in the Weissenberg geometry. In both cases, the deformation can be measured by imaging, from above, the intersection of the surface and a laser sheet projected at a low angle. aration of up to 2.5 mm between the two plates. An array of force transducers, flush mounted to the stationary plate, enabled measurement of both normal stress difference coefficients, α 1 and α 2 (Bounoua et al., 2016a). This device, by varying the plate spacing and using fibers with lengths of less than 0.5 mm, was used to characterize the normal stresses over a range of confinements from H/L = 2 to 11. Fig. 13a compares measurements at nL 2 d ≈ 3.0 of Bounoua et al. (2016a) as a function of confinement with those from the free surfaces where H/L was very large. The results are qualitatively similar, with α 2 < 0, α 1 > 0, and the magnitude of the coefficients are smaller at higher aspect ratios. The results for α 1 and α 2 at the smaller aspect ratio of A = 10 indicate that the magnitudes decline as the confinement is relaxed, approaching the values attained using the much larger geometries. For the data with A = 18, the normal stress coefficients measured with the torsional flow cell are nearly independent of the confinement. However, the magnitude of the coefficients as measured at H/L >> 1 are quite different than those for which H/L < 10.

Origin and effects of normal stresses
Normal stress differences within non-colloidal suspensions can arise from hydrodynamic interactions as well as from contact stresses. For spheres at concentrations of ϕ > 0.20, the balance of an implied shear flow, disturbance velocities created by the presence of the particles, and interparticle forces that maintain hydrodynamic interactions causes an asymmetry in the relative positions of particles pairs. This asymmetry generates a hydrodynamic contribution to the normal stress differences, and there are also contributions from the contact forces. Note that a dilute suspension of spheres (ϕ < 0.02) has no normal stress differences, as the symmetry of the velocity disturbance does not allow for it. Unlike spheres, individual fibers in a sheared fluid can induce normal stress differences depending on the instantaneous angle with respect to the flow direction. However, normal stress differences do not exist in a dilute suspension since the orientation distribution remains symmetric relative to the flow direction.
Simulations by Snook et al. (2014) demonstrate that even under concentrated conditions, symmetry of the orientation distribution about the direction of flow is not broken. Consequently, the predicted normal stresses shown in Fig. 12 originate solely from contact interactions. The predicted values of the normal stress coefficients correspond closely to the measured values. At fixed aspect ratio, the magnitudes of α 1 and α 2 are lower for low number densities since there are many fewer contacts. Likewise, fewer contacts occur at higher aspect ratios and fixed concentration, hence the normal stress coefficients (magnitudes) are lower. The relative signs of the coefficients can also be discerned from the simulations, which indicate the fibers align strongly with the direction of flow. As a result, the contact interactions, which are purely repulsive, are strongest in the gradient direction and weakest in the flow direction. This creates the positive values for α 1 and negative values for α 2 . Fig. 12 shows predictions from simulations with periodic boundaries rather than with solid bounding planes separated by a distance H in the gradient direction. However, Snook et al. (2014) also incorporated boundaries into the simulations and verified that the simulations predict that normal stress coefficients depend on confinement. This occurs because the bounding walls of the shearing geometry alter the orientation distribution, and hence the contact interactions between particles which control the normal stresses.
The existence of non-zero normal stresses in a sheared suspension of fibers implies that particles will migrate in inhomogeneous shearing flows according to theories for spherical particles (Nott and Brady, 1994;Nott et al., 2011). This migration creates a non-uniform concentration that further complicates predictions and rheological analysis of the flowing suspensions. The net migration of spheres from regions of high to low shear rates has been reported in multiple publications (Leighton and Acrivos, 1987;Altobelli et al., 1991;Hampton et al., 1997;Snook et al., Fig. 13 The normal stress coefficients as a function of confinement, H/L, for a) α 2 , the second normal stress difference relative to the shear stress, and b) α 1 , the first normal stress difference coefficient. The results of Bounoua et al. (2016a) are compared to those from the free-surface measurements Snook et al. (2014) (marked as H/L 1) at similar concentrations of nL 2 d ≈ 3.0 and aspect ratios A between 10 and 18. 2016). Currently, there are few studies of the shear-induced migration of rigid fibers, though these studies qualitatively agree with suspensions of spheres. For example, Mondy et al. (1994) measured the spatial distribution of ellipsoids of aspect ratios A = 2 to 18.4 within a wide-gap Couette flow and found results similar to that for spheres. More recent measurements (Strednak et al., 2018) indicate that fibers migrate within tube flow, away from the walls and toward the center as expected from observations of spheres. However, significant migration occurs at much smaller volume fractions in suspensions of fibers as compared to spheres. Fig. 14 shows the position dependent concentration at steady-state for two different aspect ratios for suspensions at ϕ = 0.03, whereas Snook et al. (2016) observed no particle migration for suspensions of spheres at ϕ < 0.1. Generally, Strednak et al. (2018) found that the extent of migration for the rigid fibers was similar at the same value of nL 2 d = 4Aϕ/π, regardless of the aspect ratio.

Conclusions
Various studies related to the experimental and numerical determination of the rheology have been outlined in this review. However, the reviewed studies have been limited to relatively rigid fibers of uniform size suspended under neutrally buoyant conditions in a viscous fluid. Moving away from this limited set of conditions introduces a number of possible complexities, such as Brownian motion and viscoelasticity, to the fluid mechanics and rheology of the suspension system.
Even for the relatively large size of fiber suspensions described in this review, weak and attractive colloidal forces were discussed as a possible cause of shear-thinning. For strongly colloidal suspensions, the rheology can be more complex (Mewis and Wagner, 2011). Notably for some colloidal suspensions of fibers, shear thickening rheology is observed, wherein the viscosity of the suspension increases with an increase in shear rate. Reviews by Barnes (1989) and Brown and Jaeger (2014) provide an in-depth look at the shear thickening of colloidal suspensions of both spherical and elongated particles. At this length scale, various mechanisms are claimed to cause non-Newtonian rheological behavior; more recently, fiber orientation was also attributed to causing shear-thickening (Rathee et al., 2020). Furthermore, shear thickening in colloidal suspensions of elongated particles is often related to the jamming transition (Cates et al., 1998;Liu and Nagel, 1998). However, shear thickening for colloidal suspensions should not be confused with the jamming transition for non-colloidal slurries (i.e., divergence of viscosity with increased concentration), as discussed in Section 3.2.
The stipulated condition of monodispersity of particle size in the experimental and numerical studies described in this review is far removed from real-world applications. For colloidal suspensions near the jamming transition, it is observed that at the same volume fraction, the viscosity and normal stresses measured are lower than that of monodisperse suspensions (Pednekar et al., 2018). Marti et al. (2005) performed rheological measurements for concentrated non-colloidal suspensions, with polydispersity in both shape (mixture of spheres and fibers) and particle size. Accurate and rigorous experimental measurements of the rheology in polydisperse systems continues to be lacking, with particularly few results for suspensions of fibers.
Moving past neutrally buoyant, rigid elongated particles, flexibility adds additional restorative forces to the hydrodynamics and particle-particle interactions, leading to complex non-linear, and irreversible dynamics affecting bulk rheological measurements (du Roure et al., 2019). For non-neutrally buoyant particles in shearing flows, sedimentation (or creaming) is countered only partially by resuspension of spherical particles (Acrivos et al., 1993); the resulting inhomogeneity in the spatial distribution of particles will affect the bulk measurements of rheology. For further reading, Shaqfeh (2019) collates several studies of suspensions of various non-colloidal particles in viscoelastic fluids, including rigid fibers.
Finally, there is a lack of detailed numerical and experimental work describing the microstructure at concentrations beyond the dilute limit and up to the jamming transition. The dependence of bulk rheological measurements on the changes in the microstructure of slurries of elongated particles in flow have been discussed (Franceschini et al., 2011;Snook et al., 2012;Tapia et al., 2017;Strednak et al., 2018). The changes in the microstructure are irreversible in concentrated suspensions which contradict the reversible nature of the Stokes equations, which govern their hydrodynamics. As an illustrative example to motivate more detailed measurements, results obtained from experiments performed by Tapia et al. (2017) show that for higher particle aspect ratio, volume fractions measured exceeded the limit described by numerical simulations performed by Williams and Philipse (2003), who prescribed the maximum isotropic packing fraction of spherocylinders, suggesting local or bulk ordering of the particles. Additionally, the appearance of yield stresses, a phenomenon generally associated with fiber flexibility and studied extensively in the pulp industry, was observed in slurries of rigid particles near the jamming limit, warranting deeper examination.
Even within the limits placed upon the suspension properties, the studies discussed in this review demonstrate that the rheology of slurries of rigid elongated particles is complex and frequently non-intuitive. Recent publications make it clear that care must be taken in interpreting rheological quantities for rigid fibers, especially when measured using confined geometries at large concentration. Failure to consider this latter aspect on the measurements can result in inaccurate assessments of the bulk viscosity, normal stress differences, and oscillatory rheology and dynamics. The enhancements in the measurement and modeling for these systems has been aided significantly by numerical simulations that have been validated against experiments. Continuing advances likely will rely on simulations that incorporate the additional physical mechanisms described in the previous paragraphs.