Rotational Particle Separation in Solutions: Micropolar Fluid Theory Approach

We develop a new mathematical model for rotational sedimentation of particles for steady flows of a viscoplastic granular fluid in a concentric-cylinder Couette geometry when rotation of the Couette cell inner cylinder is prescribed. We treat the suspension as a micro-polar fluid. The model is validated by comparison with known data of measurement. Within the proposed theory, we prove that sedimentation occurs due to particles’ rotation and rotational diffusion.


Introduction
The classical water-based drilling muds contain only water and clay and their performances are generally poor. Polymers used in drilling fluids improve stability and cutting removal. Currently, the polymer-based drilling fluids represent 15 to 18% of the total cost (about 1 million) of petroleum well drilling [1]. The reason is that such fluids appear to have load carrying capabilities, or, in other words, a yield stress, associated with the solid-like state and which primarily arises from the colloidal forces between the smallest suspended particles. Furthermore, in such systems, when the agitation is increased, a fluid-like state is recovered. Even in the fluid state, these materials usually show highly nonlinear behavior. The complex rheology related to the change of the internal state in the suspensions is still rather poorly understood.
Due to scale separation between colloidal and non-colloidal particles, polymer-based drilling muds with cuttings and other materials (concrete casting, foodstuff transport, etc.) can be considered as suspensions of noncolloidal particles embedded in a yield stress fluid. Substantial progress in the understanding of the behavior of such materials can thus be made by studying the impact of adding noncolloidal particles to a yield stress fluid of known properties [2]. Here, we develop a new mathematical model for suspensions embedded in a polymer-based fluid. To validate the approach, we consider rotation flows between two concentric cylinders when the external cylinder is fixed and the inner one rotates with a prescribed speed. The principal feature of such flows is the particles' migration toward the external boundary.
It is proved by Svedberg [3] that the particles' migration effect occurs also in pure colloidal suspensions. Moreover, sedimentation happens at high rotations. Currently, Svedberg's ultracentrifugation is known as an effective tool for studies of interaction between macromolecules of colloidal systems [4]. In this method, a highly disperse colloidal solution is enclosed in a wedge-shape cell rotating about an axis coinciding with the apex of the wedge [5]. Samples are centrifuged at speeds to produce sedimentation and shallow concentration gradients. A great progress in the study of such a flow was achieved by applying the diffusive Lamm equation based on the special empirical sedimentation coefficient [6]. Here, we restrict ourself to rotation flows of suspensions between two concentric cylinders paying attention to comparison with known laboratory experiments. We use methods of mechanics of a continuum by applying conservation laws only and where τ is the shear stress andγ is the shear rate. The rates' definition will be given below.
In the micro-polar fluid theory allowing for internal spins, stress tensor loses symmetry, couple stress appears, and the angular momentum equation should be included into conservation laws. Instead of Equation (1), one writes for simple shear flows the following equation: whereγ a is a shear rate related to particle rotation, and η sk is an additional viscosity. In what follows, we call η sk skew-symmetrical viscosity motivated by rheological relationships which will be discussed later.
When the volume concentration of particles φ is equal to zero, the viscosity η s of the interstitial fluid is assumed to be known. Clearly, η s is the usual shear viscosity and it is well known how it can be measured. The skew-symmetric viscosity manifests itself when φ > 0, and poor knowledge of it hampers applications of the micro-polar fluid models. The above description of possible applications supports the view that it is of importance to know how skew-symmetric viscosity depends on the particle concentration. The goal of the first sections of the present paper is to get a better insight of how η sk depends on φ.
In Sections 1-5, we validate some formula for the skew-symmetric viscosity and apply it in Section 6 to the problem of particles migration for the steady flows of yield stress granular fluids in a concentric-cylinder Couette geometry when rotation of the Couette cell inner cylinder is prescribed.

Micro-Polar Fluids
Here, we remind notions concerning a micropolar fluid within the theory of the Cosserat continuum. Such a fluid exhibits micro-rotational effects and micro-rotational inertia; the fluid can support couple stresses and body couples and possesses a non-symmetric stress tensor. The theory of micro-polar fluids goes back to [23][24][25], where gyration tensor, inertial spin, conservation of micro-inertia, and objectivity of micro-deformation rate tensors are derived and discussed. For an overview of developed theories, we refer the reader to [26,27].
In the Cosserat continuum [28], each material point is treated as a rigid body in the following sense. To such a point with the Lagrange coordinates (t, ξ), one can assign a position vector x(t, ξ) in the three-dimensional Euclidean space and three orthogonal directors d i (t, ξ), i = 1, 2, 3. Rotation of the vectors d i is governed by a rotation orthogonal tensor Q(t, ξ). The rotation velocity tensor is skew-symmetric, and it enjoys the representation formula where ω(t, x) is the micro-rotational velocity vector, 2ω = e i × (Ω · e i ) = : Ω.
Here, {e i } 3 1 is any orthogonal basis in R 3 , and is the Levi-Civita third order tensor, a, b, c = a · (b × c), e i × e j = sij e s , sij ≡ e s , e i , e j , ( : Given a 3 × 3-matrixes A and B, we use the notation A * for the adjoint matrix such that and the scalar product A : B is defined by A : With v(t, x) standing for the velocity of the mass center of the Cosserat material point (t, ξ), the micropolar fluid is characterized by two rates of strain tensors B and A: Here, we use the notations (∇v) ij = ∂v i /∂x j , (∇v) * ij = ∂v j /∂x i . An instant stress state of such a fluid is characterized by the couple stress tensor N(t, x) in addition to the Cauchy stress tensor T(t, x). Let S stand for the viscous part of the stress tensor, T = −p I + S. In what follows, we use the symmetric and skew-symmetric parts of a matrix D: In the typical stress-strain relation of Newtonian or non-Newtonian fluids, the scalar factor η s is the viscosity, p is the pressure, with I and T being the identity and stress tensors. The tensor T is symmetric as is well known in the classical fluid mechanics theory, i.e., T * = T. Such a symmetry results from momentum and angular momentum laws. On the other hand, the angular momentum law is valid automatically if we postulate symmetry of T. This is why nobody invokes such a law in applications. We remind readers that the constitutive laws of a simple micropolar fluid are [26] where η s , η sk are the viscosities and γ is angular viscosity. The first rheological equation in (5) suggests that the contributions of the symmetric part B s = (∇v) s and skew-symmetric part B a of the rate of strain tensor B into local stress state are different. It is proved in [25] that both the rate of strain tensors B and A are objective. Let us introduce the relative angular velocity Observe that in fluid mechanics shear stress and shear rates in (1) are defined as follows: Thus, the skew-symmetric viscosity characterizes how relative micro-rotations contribute into the local stress state.
Observe that in the Cosserat continua the Cauchy stress tensor is not symmetric and the vector t = e i × (T · e i ) = : T is a stress symmetry defect measure in the sense that the equality t = 0 implies T * = T and vice versa. By definition of t, we have the formula The momentum and the angular momentum laws are where ρ is the density, J is the micro-inertia scalar, f is the mass force vector and The density ρ satisfies the mass conservation laẇ

Skew-Symmetrical Viscosity of Dilute Suspensions of Rigid Particles
Let us consider Couette-like steady flows of suspensions between two parallel planes in the x-direction when the upper plane y = h is fixed and the lower plane y = 0 moves in the x-direction with the velocity V. The volume particle concentration φ is assumed to be fixed.
We outline the method for determination of the skew-symmetrical viscosity η sk . Assume that the stress applied to the moving plate can be measured. By continuity, one can tell the fluid stress in the nearby fluid region. On the other hand, one can calculate such a fluid stress by one or another mathematical model. First, we determine the stress S m 1 at the moving plane by applying the micro-polar fluid theory. Then, we calculate the stress S v 1 at the moving plane by applying the Navier-Stokes theory. We equate the two stresses and derive the skew-symmetrical viscosity η sk from the equality S m 1 = S v 1 . Let us first treat the suspension as a micropolar fluid with the prescribed volume particle concentration φ. The above assumptions upon the flows suggest that the velocity vector v, the micro-rotational velocity vector ω and the pressure p depend on the vertical variable y only: For such flows, the matrices ∇v and Ω become ∇v = Hence, Projections of the momentum and the angular momentum equations on the x and z-directions become respectively, where the tensor components are given by the formulas Constitutive laws (5) take the form We note that system (10)-(11) does not contain pressure. As is well known, it can be restored from the momentum equation projected on the y-direction.
The viscosity η sk (φ) vanishes when φ → 0. The same is true for the relative viscosity which we represent via the expansion series Whereas the velocity v satisfies the no-slip boundary conditions, we require that where α(φ) = α 0 φ. The latter condition implies that the micro-rotations agree with macrorotations at the boundary [29]. For the Couette-like flows, we arrive at the following boundary-value problem in the domain 0 < y < h: 2γ We solve the boundary value problem (14)-(16) looking for (v, w) as the asymptotic expansion series Setting these series in (14)- (16), one can write each equality in (14)- (16) in the form The coefficients v i and ω j are determined from the conditions (· · · ) k = 0 for any k.
Particularly, if k = 0, we derive the following boundary value problems for the functions v 0 (y), ω 0 (y): Similarly, we find that the function v 1 satisfies the boundary value problem ∂ ∂y Solving these problems, we find that Starting from the definitions (12) related to the micro-polar fluid theory, we can write the expansion series for the stress S 12 (φ) as follows: Clearly, Now, we can calculate the relative stress at the moving plane: It is assumed that both the stresses S m 12 (φ) and S m (0) 12 are measured at the same velocity V of the moving plate.
Let us consider flows within the same Couette geometry, starting from the Navier-Stokes theory. In such a theory, the stress tensor S is symmetric. Denoting S = S 12 , one can find the velocity v(y) by solving the boundary-value problem where η(φ) is the effective viscosity. Let S v 1 stand for the stress S at the moving plate, S v 1 = S| y=0 . Given S v 1 , one can derive from (22) the following formula for the apparent viscosity η(φ): Observe that S v 1 has negative values since v(y) is a decreasing function of y. It follows from (23) It is assumed that both the stresses S v 1 (φ) and S 1 (0) v are measured at the same velocity V of the moving plate.
For dilute suspensions, the left-hand side of (24) is given by the Einstein formula [30] We equate the relative stresses: (24) and (25) that Hence, Λ = E. By the above arguments, we conclude that the skew-symmetric viscosity for dilute suspensions satisfies the representation formula where E 2.5 is the Einstein factor.

Remark 1.
By the same asymptotic arguments, we can conclude that for the Couette flows between two concentric cylinders formula (26) becomes where G is the geometrical factor equal to (R 2 /R 1 ) 2 , with R 1 being the smaller radius.
One more conclusion from the above arguments is that there is a correlation between the Navier-Stokes theory and the micro-polar fluid theory: where η s , η sk are the micro-polar fluid viscosities of the suspension and η(φ) is the apparent viscosity of the same suspension described by the Navier-Stokes rheology. The law (28) is verified by the asymptotic series argument for dilute suspensions, with the left-hand side given by the Einstein law On the other hand, there is an extended Krieger-Douhgerty empirical closure [31] for dense suspensions, where φ * is a maximal volume concentration. Such a closure suggests that, by setting (29) in (28), we can define the skew-symmetric viscosity η sk (φ) as follows: In the next sections, we verify this empirical formula by studying flows of dense suspensions paying attention to particle rotation.

Flows of Suspensions of Rigid Particles in the Herschel-Bulkley Fluid
By the theoretical approach of Chateau et al. [32], Ovarlez et al. [9] showed that the flows of yield stress suspensions in a concentric-cylinder Couette geometry can be modeled by a Herschel-Bulkley behavior of same index as their interstitial fluid. The theory was proved to be in an agreement with the laboratory experiments based on the magnetic resonance imaging techniques [9].
We are going to address the same experiments as in [9] to verify formula (30). First, we extend the constitutive laws (5) to allow for the yield stress rheology. According to [33], the Cosserat-Bingham fluid rheology is defined as follows: where and τ * and τ n are yield stresses; the unknown plug tensors S p and N p obey the restrictions It is proved in [34] that the formulation (31) Similarly, the constitutive law (32) is equivalent to the inclusion N ∈ ∂V n (A) with V n (D) = γ|D| 2 + τ n |D|, ∀D ∈ R 3×3 . The meaning of the plug zone |N(x, t)| ≤ τ n is discussed in [35].
Let T 0 be a characteristic time. We denote the dimensionless second invariant of the rate of strain tensor B 0 by I: I = T 0 |B s |. To transform the Cosserat-Bingham fluid constitutive laws (31) and (32) into the Cosserat-Herschel-Bulkley fluid rheological equations, we assume that Observe that the classical Hershel-Bulkley model results from the constitutive laws (31)-(33) if the particle volume concentration φ vanishes.
We consider steady axially symmetric flows of a suspension between two coaxial cylinders centered on the z-axis. The inner cylinder of the radius R 1 rotates with the angular velocity Ω 0 [s −1 ] and the external cylinder of the radius R 2 is fixed. The volume particle concentration φ is assumed invariable along the radial coordinate. The case of variable φ will be addressed in Section 6.
In what follows, we use the unit vectors e r , e ϕ , e z of the cylindrical coordinate system. The assumption of axially symmetry of flows suggests that the velocity vector v, the microrotational velocity vector ω and the pressure p depend on the radial variable r only: v = v(r)e ϕ , ω = ω(r)e z , p = p(r), rot v = ∂v ∂r For such flows, the matrices ∇v and Ω in the cylindrical coordinate system become ∇v = where we denoted ∂v/∂r by ∂v for simplicity. Here, Hence, As for the rate of strain tensor A, we have that A zr = ∂ω/∂r and A ij = 0, otherwise. Projections of the momentum Equation (7) respectively, where the tensor components are given by the formulas S ϕr = e ϕ · S e r , S rϕ = e r · S e ϕ , N zr = e z · N e r .
Given a vector e, we apply the notation (S e ) i = S ij e j . Observe that the other components of S and N are equal to zero. In what follows, we use the equation resulting from projection of the momentum equation (7) onto the vector e r . We calculate that The constitutive laws (31)-(32) become The boundary condition (13) for the angular velocity ω and the no-slip condition for v become To study numerically the boundary-value problem (33)-(42) in the annulus R 1 < r < R 2 , we pass to dimensionless variables: Observe that the dimensionless yield stress τ * 1 is the inverse of the Bingham number for the Couette flows: In new variables,

Skew-Symmetric Viscosity versus Particles Concentration
Here, we apply the mathematical model developed in the previous section to justify formula (30) for the skew-symmetric viscosity. To perform calculations, we fix parameters of the interstitial Hershel-Bilkley fluid. Rheological constitutive law of such a fluid results from (31) by setting φ = 0: We denote τ y = The concentrated emulsion obeying Equation (53) was considered in [9] with τ y = 22 [Pa] and η HB =5.3 [Pa·s 1/2 ]. Given the angular velocity Ω [rpm] of the rotating inner cylinder, we calculate that T 0 = 1/Ω 0 = (60/Ω)[s]. We consider the same fluid as in [9], hence one can define the consistency η s0 and the yield stress τ * as follows: Observe that consistency depends on Ω because of the special choice the characteristic time T 0 and the definition of the dimensionless invariant I of the rate of strain tensor B s .
It is proved in [35] that the rotation yield stress τ n causes the appearance of clusters of particles, with each cluster being a plug zone which rotates as a rigid body. Conglomerates of particles were not observed in [9] for the Couette flows of suspensions between two rotating cylinders; this is why τ n and τ n1 can be neglected. As for the dimensionless angular viscosity γ 1 and the boundary-value dimensionless parameter α 0 , we variate them to fit experimental data.
Approximate solutions of the system (43)-(50) can be obtained by regularization [36]. Given a small positive δ, we substitute the dimensionless invariant I in (44) and (45) by I δ , where with δ 0. First, we tune the model (43)-(50) by setting φ = 0 and addressing the pure interstitial Herschel-Bulkley fluid with τ y (0) = 22[Pa] and η HB (0) = 5.3[Pa · s 1/2 ] as in [9]. Equations become Observe that in such a case the model (55)-(59) depends on one parameter τ * 1 (0) only. Given the angular velocity Ω[s −1 ], we find from (52) the value of the Bingham number Bn(0) by the formula Figures 1 and 2 depict very good agreement of calculation results for φ = 0 with experiment data [9] for different angular velocities Ω 0 [s −1 ] but in the case of passage to the effective Bingham number Bn e (0): We think that such a discrepancy between the measured and effective Bingham numbers is due to the following reasons. Real 3D-flows are described by 1D-equations. The gravitation, the height of the annulus region, and the lateral boundaries effect are not taken into account. It may be that viscoelastic fluid property is also of importance, which calls for more adequate modeling. However, the dimensionless yield stress τ * 1 (0) stipulates the size of the plug zone and, at the same time, it varies inversely with the angular velocity Ω. Now, we consider suspensions assuming that they can be modeled by a Herschel-Bulkley behavior of the same index as their interstitial fluid, with their consistency and their yield stress depending on the particle volume fraction. To determine the function τ * 1 (φ) in Equations (55)-(59), we apply correlations proposed in [9] for n = 1/2: In our notations (see (33)), equality (63) becomes It follows from the definitions (60) and (61) that With the function τ * 1 (φ) given by (64), we solve Equations (55)-(59) and find an agreement with experiment data [9]. Calculations and laboratory data depicted in Figure 3 are related to φ = 0.3 for Ω = 2, 5, 10, 20, 50 and 100 [rpm].
Let us comment on some discrepancy between calculations and data of measurement. One can see in Figures 4b and 7 that, with increasing the rotational velocity at a constant volume concentration φ = 0.1, the theoretical results are shown to better agree with experiment. It is a problem of calculations. The reason is that we consider the viscoplastic fluid and the plug zone (with zero velocity and low shear stress) near the external cylinder being larger at low rotations. The governing system of equations becomes degenerate in such a zone. Mathematical theory of degenerate systems of differential equations and corresponding numerical methods built into Wolfram Mathematica are not well developed yet. There is one more difficulty related to the Hershel-Bulkley fluid viscosity (with the power n = 1/2) becoming infinite when the velocity gradient vanishes somewhere. To get over these difficulties, we apply the regularization technique and substitute the invariant I of the rate of strain tensor B 0 in Equations (38)-(40) by its non-vanishing approximation I δ .
A comparison of the results in Figure 1 (upper curve), and Figure 4 at the rotational velocity Ω = 100 rpm shows that, with increasing the volume concentration, the calculated results better agree with experiments at intermediate concentrations φ=0.1; at lower (φ = 0) and at higher concentrations (φ = 0.3), the deviations increase. The reason is that there are no data of measurement for the dimensionless parameters γ 1 and α 0 . To choose them, we apply the method of least squares based on minimizing the functional F(γ 1 , α 0 ). As it happened, a discrepancy between the measured data and calculations for different concentrations and angular velocities is due to the optimal choice of these unknown parameters.
The above calculations confirm that Equation (30) for the skew-symmetric viscosity η sk (φ) can be of use.
Observe that comparison with experiments for colloidal suspensions is contained in Figures 1, 2 and 5 since, in the case φ = 0, the suspension becomes a pure colloidal fluid. Figure 1. Calculated dimensionless velocity profiles (solid lines) versus the radial variable for pure interstitial Herschel-Bulkley fluid without particles, φ = 0. The lines from the bottom upward correspond to Ω = 2, Ω = 5 and Ω = 100, respectively. Stars, balls, and triangles stand for measurement data [9] in the cases Ω = 2, Ω = 5 and Ω = 100, respectively.

Rotational Sedimentation
Here, we consider more general mathematical model allowing for non-uniform particle distribution. We introduce the mass concentration of particles as follows: where ρ is the total density,ρ s is the particle density andρ f is the density of the interstitial fluid. Given c, one can restore from (65) the volume concentration and the total density by the formulas Due to these formulas, any given function of the volume concentration like the relative viscosity ε(φ) can be defined in terms of the mass concentration c. It is explained in [7] that c satisfies the conservation law where l is the concentration flux obeying the generalized Fick equation The scalar parameters D[cm 2 /s], D p [cm 3 · s/g], and D ω [cm 2 · s], stand for the diffusion, barodiffusion, and spin diffusion coefficients.
Observe that, instead of (5), the couple stress tensor N is prescribed by the rheological equation to meet the entropy production law [7] where the skew-symmetric matrix : a is defined by the formula ( : a) ij = a k ikj , a = a i e i , in any orthogonal basis {e i } 3 1 . It is due to spin diffusion that the Ségre-Silberberg effect is explained within the micropolar theory [7].
Due to the identity rot ω × b = 2(∇ω) a · b, ∀b ∈ R 3 , it follows from (68) and (69) that in agreement with the definition of rotary diffusion [37]: "Just as the translational diffusion coefficient is calculated in terms of the drag force, so the rotary diffusion coefficient is expressed in terms of the moment of the forces on a particle executing a rotary movement." It follows from (9) that for steady flows the mass conservation law becomes divv = 0.
With the above definitions, we arrive at the following conservation laws for steady flows: Under the assumption that c = c(r), we arrive at the formula ∇c = e r ∂c/∂r. Due to equation divv = 0, we obtain that, for the rotation flows, the equation div(ρcv) = 0 holds. Now, it follows from (72) that 0 = div l = 1 r ∂(rl) ∂r and rl = const.
At the same time, the now-flow boundary conditions l · n| R i = 0 imply that l = 0 and l = 0. Thus, the particles mass concentration obeys the equation l = 0 or ∂c ∂r Due to (37), the latter equation can be written as In what follows, we use the representations D p = D * p c(1 − c), D ω = D * ω c(1 − c) since both D p and D ω vanish at c = 0 and c = 1. Given a mean value c 0 of c, we set the following condition: As for the particle migration, there is one more approach based on the Fick law. This is known as the Lamm equation for a highly disperse colloidal solution enclosed in a wedge-shape cell rotating at an angular velocity Ω about an axis coinciding with the apex of the wedge [5].
Let us return to Equation (73). One more consequence of the equality l = 0 is that the rheological Equation (69) reduces to N = 2γA in the Couette geometry.
We introduce dimensionless parameters In new notations, Let us pass to dimensionless variables. Then, Equations (73) and (74) become We summarize the governing equations as follows. We look for the dimensionless functions v (r ), ω (r ) and c(r ) which satisfy Equations (43)-(50), (75) and (76). Observe that these equations are not decoupled since the relative viscosity ε(c) in Equations (44) and (45) depends on the particle concentration c. Figure 11 depicts results of calculations of concentration along the radial variable. We apply the Wolfram Mathematica solver for ordinary differential equations. Agreement between calculations and experiment [9] is achieved for Ω = 10 2 rpm with the choice D * ω /(2D) = 5 × 10 −5 [s 2 ]. The sedimentation effect happens when we increase Ω to 14 × 10 3 [rpm]. We prove that such an effect is due to the rotational diffusion D ω since the particle separation does not happen when D ω = 0.

Remark 2.
It is known that, for many practical purposes, the polymers or colloidal particles can be regarded as rigid particles. Examples are triblock Janus particles which can be modeled as crosslinked polystyrene spheres whose poles are patched with sticky alkyl groups, and their middle band is covered with negative charges [38]. This is why the above results on rotational sedimentation can also be applied to polymer flows. As was proved by Svedberg, such flows are of great importance in the studies of the polymer structure. Our contribution is that we propose an alternative approach to the Lamm equation based on the empirical notion of the sedimentation coefficient [6]. The advantages are that we apply the conservation laws of continuum mechanics and take into account the shape of the particles. Moreover, we prove that the polymer sedimentation is due to its rotation. This result suggests a new direction of laboratory studies on polymer flows. Figure 11. Profiles of mass concentration c. Both the solid line that is based on calculations and dots standing for experimental data [9] correspond to Ω = 10 2 rpm. Agreement between calculations and experiment is achieved by the choice D * ω /(2D) = 5 × 10 −5 [s 2 ]. The dashed line corresponding to calculations reveals the sedimentation effect when we increase Ω to 14 × 10 3 [rpm].

Discussion
We address the rotational sedimentation of particles for steady flows of yield stress granular fluids in a concentric-cylinder Couette geometry. Apart from the Lamm equation approach, we do not use the empirical sedimentation coefficient. Instead, we apply conservation laws of the micro-polar equations which allow for particle rotation. We prove that it is due to the rotational diffusion that the particle sedimentation occurs at high angular velocity of the Couette cell inner cylinder. To validate the mathematical model, we perform a comparison with published data of measurements by choosing the relative viscosity related with the particle rotation. First, we justify analytically this choice for dilute suspensions starting from the Einstein correlation for the apparent viscosity. As for dense suspensions, we apply the Krieger-Douhgerty empirical closure for the apparent viscosity. Though we performed calculations for steady flows, the developed approach allows for unsteady flows and non-spherical particles due to the micro-inertia tensor involved into the angular momentum conservation law.