J Meibohm

F Candelier

T Rosén

J Einarsson

F Lundell

B Mehlig

Angular velocity of a spheroid log rolling in a simple shear at small Reynolds number

Keywords: numbers: 83.10.Pp, 47.15.G-, 47.55.Kf, 47.10.-g

We analyse the angular velocity of a small neutrally buoyant spheroid log rolling in a simple shear. When the effect of fluid inertia is negligible the angular velocity ω equals half the fluid vorticity. We compute by singular perturbation theory how weak fluid inertia reduces the angular velocity in an unbounded shear, and how this reduction depends upon the shape of the spheroid (on its aspect ratio). In addition we determine the angular velocity by direct numerical simulations. The results are in excellent agreement with the theory at small but not too small values of the shear Reynolds number, for all aspect ratios considered. For the special case of a sphere we find

where s is the shear rate and Res is the shear Reynolds number. This result differs from that derived by Lin et al. [J. Fluid Mech. 44 (1970) 1] who obtained a numerical coefficient roughly three times larger.

I. INTRODUCTION

The angular dynamics of a small, neutrally buoyant particle in a viscous fluid is determined by the local fluidvelocity gradients. In the creeping-flow limit the particle rotates in such a way that the torque exerted by the fluid vanishes at every instant in time. For larger particles this is no longer true, in general. The particle angular momentum may change as the particle rotates, and the fluid near the particle is accelerated. It is usually a difficult problem to compute the effect of particle and fluid inertia on the dynamics of particles in flows.

In this paper we calculate the leading inertial correction to the angular velocity of a neutrally buoyant spheroid that rotates in an unbounded shear so that its axis of symmetry is aligned with the flow vorticity (Fig. 1). This motion is referred to as 'log rolling'. For an oblate spheroid this orbit is stable, whereas it is unstable for a prolate spheroid [START_REF] Einarsson | Rotation of a spheroid in a simple shear at small Reynolds number[END_REF]. The calculations summarised in this paper provide the corresponding angular velocity.

In the log-rolling orbit the locus of the particle surface does not change as a function of time. This simplifies the calculations considerably, because we can consider the steady problem. The latter is also easiest to analyse in direct numerical simulation because it is not necessary to change the numerical mesh as the particle moves.

In the log-rolling orbit the angular velocity aligns with the vorticity axis ê3 (Fig. 1), and

ω ≡ ω • ê3 ∼ - s 2 + 0.0540 3sD 10π Re 3/2 s ( 1
)
to order O(Re 3/2 s ). The first term on the r.h.s. corresponds to the angular velocity in the creeping-flow limit, equal to -s/2 for an unbounded shear with shear rate s [START_REF] Jeffery | The motion of ellipsoidal particles immersed in a viscous fluid[END_REF]. The second term on the r.h.s describes the inertial correction for small shear Reynolds number Re s ≡ a 2 s/ν, and ν the kinematic viscosity of the fluid. Furthermore D is a shape parameter that depends on the aspect ratio Spheroid log rolling in a simple shear. The r2coordinate is aligned with the shear direction, and the r3coordinate is the negative vorticity direction. The axis of symmetry of the particle and its angular velocity ω = ωê3 are aligned with the flow vorticity.

of the spheroid. For a sphere we have D = 10π/3, and the general expression is given in Eq. (20) below. It turns out that D > 0 for all aspect ratios, so that the effect of fluid inertia reduces the magnitude of the angular velocity. But we note that there are other examples where fluid inertia increases the angular velocity of a neutrally buoyant sphere [START_REF] Bluemink | A sphere in a uniformly rotating or shearing flow[END_REF].

For a spherical particle Eq. ( 1) reduces to ω/s ∼ -1/2 + 0.0540 Re 3/2 s . This is different from a result obtained by Lin et al. [START_REF] Lin | Simple shear flow around a rigid sphere: inertial effects and suspension rheology[END_REF] who found a coefficient of 0.1538 for the inertial correction, instead of 0.0540. Our result is consistent with that of Stone, Lovalenti & Brady [START_REF] Stone | Inertial effects on the rheology of suspensions and on the motion of individual particles[END_REF] who used a different technique. We also performed direct numerical simulations of the problem to determine the angular velocity. For a spherical particle the results are shown in Fig. 2. We see that Eq. (1) agrees well with the simulation results when Re s is neither too large nor too small.

Why do the numerical results agree with Eq. ( 1) only at intermediate values of Re s ? This is a consequence of the fact that fluid inertia causes a singular perturbation of the creeping-flow problem. Even at very small shear Reynolds number the perturbation is not negligi- Re s ! ! (0) ()

! (0) ()+ 1 2  FIG. 2. Re 3/2
s -correction to the angular velocity of a spherical particle. Open symbols show results of our direct numerical simulations (Section V) for κ = 0.01 ( ), 0.025 (♦), and 0.05 (•). Here κ = 2a/L, a is the radius of the sphere, and L is the linear dimension of the simulation domain. The thick solid (red) line shows Eq. ( 1), valid for an unbounded shear. The thin solid (blue) lines correspond to a quadratic Resdependence, they show fits to Eq. ( 49). Also shown are experimental results by Poe & Acrivos (Section VI), filled triangles. The creeping-flow limit of the angular velocity in the finite, bounded system is denoted by ω (0) (κ), and ω (0) (κ) → - 1 2 as κ → 0. The inset shows ω (0) (κ)+ 1 2 as a function of κ ( ). The thick solid line (red) shows Eq. ( 48), with the fit-parameter C = 0.22. ble far away from the particle, at distances larger than the Saffman length S ≡ a/ √ Re s . This means that a regular expansion of the solution to the particle (the 'inner solution') must match an approximate 'outer solution' that takes into account convective fluid inertia but does not fulfill the no-slip boundary conditions on the particle surface. The power Re 3/2 s in Eq. ( 1) is due to the singular nature of the perturbation. The small parameter of the problem is ≡ √ Re s .

The numerical simulations shown in Fig. 2 were performed for a bounded shear, the simulation domain is a cube with side length L (and we define κ ≡ 2a/L). We expect Eq. ( 1) to agree with the numerical results when L S , in the opposite limit the perturbation is in effect regular. Fig. 2 confirms this picture: the numerical results converge to Eq. (1) as L increases, but there are substantial deviations at small values of Re s , as mentioned above. In this regime the inertial correction to the angular velocity is quadratic in Re s (thin solid blue lines in Fig. 2). That there is no term linear in Re s is a consequence of the symmetry of the problem. Naturally the perturbation theory must also fail at large Re s . Eq. ( 1) works reasonably well up to Re s ≈ 5 × 10 -2 .

For a spheroidal particle the steady-state angular velocity depends on the particle shape. Fig. 3 shows our s -correction to the angular velocity of a neutrally buoyant spheroid log rolling in a simple shear, as a function of the aspect ratio λ. Shown are results of our direct numerical simulations (Section V, ) for Res = 3 × 10 -3 and κ = 0.01. The theoretical result [Eq. ( 1), valid for an unbounded shear] is shown as a solid red line. The red arrow points to the corresponding axis. Also shown is the surface area |S | of the particle (black dashed line). The corresponding axis is on the r.h.s. (black arrow). results for the inertial correction of the angular velocity as a function of the particle aspect ratio at small but finite Re s . For a prolate spheroid the aspect ratio is defined as λ ≡ a/b, where a is the major semi-axis length of the particle, and b is the minor semi-axis length. For oblate spheroids the aspect ratio is defined as λ ≡ b/a [START_REF] Einarsson | Rotation of a spheroid in a simple shear at small Reynolds number[END_REF]. Since the shear Reynolds number is defined in terms of the major semi-axis length of the particle we leave a unchanged as we vary the particle shape.

In Fig. 3, the theoretical result (1) (red line) is compared to results of direct numerical simulations (white squares) of a large system (κ = 0.01) at Re s = 3 × 10 -3 . We find excellent agreement.

The asymptotic matching method described in this paper is slightly different from that used by other authors, and has distinct advantages for the present problem. Therefore, we briefly comment on the differences here. We obtain the outer solution as an expansion in kspace, applying the approach described in Ref. [START_REF] Candelier | Note on the method of matched-asymptotic expansions for determining the force acting on a particle[END_REF]. But we use a new way of generating the terms in this expansion that substantially simplifies the asymptotic matching. We obtain the angular velocity to order 3 in terms of a lower-order solution, valid to order . The method employed by Stone, Lovalenti & Brady [START_REF] Stone | Inertial effects on the rheology of suspensions and on the motion of individual particles[END_REF], based on the reciprocal theorem, has the same advantage. Our calculations show that their method yields precisely the same integral expressions as our asymptotic matching. It remains to be investigated how general this correspondence is.

II. FORMULATION OF THE PROBLEM

We consider a small neutrally buoyant spheroid rotating in a simple shear with shear rate s. We assume that the axis of symmetry of the particle is aligned with the fluid vorticity (Fig. 1) and that the centre-of-mass of the particle is advected by the flow. Our aim is to determine how weak inertia affects the angular dynamics. In the log-rolling orbit angular velocity ω aligns with the vorticity axis ê3 at all times. We compute ω ≡ ω • ê3 as an expansion in ≡ √ Re s :

ω = ω (0) + ω (1) + 2 ω (2) + 3 ω (3) + . . . . (2)
We calculate ω from the condition that the torque must vanish in the steady state, τ = 0. The torque is determined by the solution of Navier-Stokes equations:

written in dimensionless variables. As time scale we take the inverse shear rate s -1 , as length scale we take the length a of the major axis of the particle, and as velocity and pressure scales we take as and µs, respectively. Here µ is the dynamic viscosity of the fluid. The torque that the fluid exerts upon the particle is measured in units of µsa 3 . All equations in the remainder of this article are written in dimensionless variables. In Eq. ( 3), u i are the components of the fluid velocity and p is the pressure. We use the implicit summation convention that repeated indices are summed over from 1 to 3. The boundary conditions are no slip on the surface S of the particle, and that the fluid velocity is undisturbed infinitely far away:

u i = ε ijk ω j r k for r ∈ S , and u i → u (∞) i as r → ∞. (4)
Here u (∞) i are the components of the undisturbed fluid velocity, r = |r|, r is the Cartesian coordinate vector with components (r 1 , r 2 , r 3 ), and the origin of the coordinate system is located at the centre of the particle. Further ε ijk is the Levi-Civita symbol. For the simple shear shown in Fig. 1 we have

u (∞) i = A (∞) ij r j with A (∞) ij = δ i1 δ j2 . (5)
Here

A (∞) ij ≡ ∂ j u (∞) i
are the undisturbed fluid-velocity gradients, and δ ij is the Kronecker symbol.

We decompose the fluid velocity into the undisturbed velocity u (∞) i and the disturbance velocity u i , and decompose the pressure in a similar way:

u i = u (∞) i + u i , p = p (∞) + p . (6)
Since the locus of the particle surface does not change as a function of time way may consider the steady disturbance problem:

2 (A (∞) jk r k ∂ j u i +A (∞) ij u j +u j ∂ j u i ) = -∂ i p +∂ j ∂ j u i . (7)
The boundary conditions for u i follow from Eq. ( 4):

u i = ε ijk ω j r k -u (∞) i for r ∈ S , (8a)
u i → 0 as r → ∞ . (8b)
Setting = 0 in Eq. ( 7) corresponds to the Stokes limit. The disturbance flow in this limit, u

i , is well-known. At finite values of , by contrast, Eq. ( 7) is difficult to solve because the inertial terms on the l.h.s. of Eq. ( 7) constitute a singular perturbation of the Stokes problem.

The singular nature of the perturbation is verified as follows. The leading large-r term in the Stokes solution of our problem decays as u (0) i ∼ r -2 . This allows us to estimate the terms on the l.h.s. of Eq. ( 7) as r → ∞:

2 A (∞) jk r k ∂ j u (0) i ∼ 2 A (∞) ij u (0) j ∼ 2 r -2 , ( 9a
) 2 u (0) j ∂ j u (0) i ∼ 2 r -5 . ( 9b
)
The inertial terms (9a) balance the viscous term

∂ j ∂ j u (0) i ∼ r -4 (10)
at the Saffman length S (equal to 1/ in dimensionless variables). Thus even at very small shear Reynolds number the effect of convective fluid inertia cannot be neglected at distances larger than the Saffman length. We conclude that the problem is singular, perturbation theory requires matched asymptotic expansions [START_REF] Proudman | Expansions at small Reynolds numbers for the flow past a sphere and circular cylinder[END_REF][START_REF] Saffman | The lift on a small sphere in a slow shear flow[END_REF]. In the 'inner region' (r < S ) regular perturbation theory can be used. The inner solution satisfies the boundary condition (8a) close to the particle. At r ∼ S the inner solution is matched to a solution obtained in the 'outer region' (r > S ) that takes into account convective fluid inertia and satisfies the boundary condition (8b).

III. METHOD A. Inner solution

We seek regular perturbation expansions of the disturbance velocity and the pressure in the inner region:

u in,i (r) = u (0) in,i (r) + u (1) in,i (r) + 2 u (2) in,i (r) + . . . , (11a) p in (r) = p (0) in (r) + p (1) in (r) + 2 p (2) in (r) + . . . . (11b)
The boundary conditions (8a) read order by order:

u (0) in,i (r) = ε ijk ω (∞) j r k -u (∞) i for r ∈ S , (12a) u (n) in,i (r) = ε ijk ω (n) j r k for n ≥ 1 and r ∈ S . (12b)
These boundary conditions do not suffice to determine the inner solution. Therefore each term in the expansion (11a) is matched to its counterpart in the outer solution at distance r ∼ S from the particle.

B. Outer solution

To determine the outer solution we solve

2 (A (∞) jk r k ∂ j u o,i +A (∞) ij u o,j ) = -∂ i p + ∂ j ∂ j u o,i +f i (13)
in the outer region. The subscript 'o' stands for 'outer'. In Eq. ( 13) we have neglected the non-linear convective term u j ∂ j u i . At the end of this Section we verify that the error is of order O( 4 ). The source term f i is introduced to account for the presence of the particle:

f i (r) = D ij ∂ j δ(r).
Here δ(r) is the three-dimensional Dirac delta function. The coefficients D ij must be determined so that the outer solution matches the inner solution. We expect that the disturbance flow due to the particle can be represented in the outer region by a (symmetric) stresslet in the ê1 -ê 2 plane. Hence we make the ansatz D 11 = D 22 = D 3i = 0 and D 21 = D 12 ≡ D (the 'dipole strength' D is determined by matching below).

Since Eq. ( 13) is linear it can be solved by Fourier transform. In Section IV we compute the expansion of the Fourier transform û o,i (k) in powers of ≡ √ Re s :

û o,i = Ti (0) + Ti (1)
+ 2 Ti

(2)

+ 3 Ti (3) + . . . . (14)
It turns out that the odd terms in this expansion are generalised functions of k [START_REF] Candelier | Note on the method of matched-asymptotic expansions for determining the force acting on a particle[END_REF]. Transforming back to real space we find:

T i (0) (r) = - 15D 20π r 1 r 2 r i r 5 , (15a)
T i (1) (r) = 0 , (15b)
T i (2) (r) = - 15D 240π 1 r 5 r 2 1 r 2 2 r i -r 1 r 2 r 2 ε ij3 r j (15c) - 1 r -r 1 δ i1 1- r 2 1 3r 2 +δ i2 r 2 1+ r 2 2 3r 2 -δ i3 r 3 3 3r 2 , T i (3) (r) = 3D 10π A ij r j . (15d)
Details are given in Section IV. Eqs. ( 15) can be used to verify that the non-linear convective term u j ∂ j u i can be neglected in the outer problem (13) to order 3 . The argument goes as follows. Eq. ( 15) shows that T i (n) (r) scale as ∼ r n-2 for n = 0, . . . , 3 in the matching region. In this region, r is of order of -1 . Using ( 14) it follows that all terms in Eq. ( 13) are of order O( 4 ) in the matching region:

2 A (∞) jk r k ∂ j T i (n) +A (∞) ij T j (n) ∼ ∂ j ∂ j T i (n) ∼ O( 4 ) (16)
for all n = 0, . . . , 3. Using the same arguments we can also estimate the magnitude of the neglected non-linear convective term in the matching region:

2 T j (n) ∂ j T i (n) ∼ O( 7 ) . (17)
It turns out that the non-linear convective term remains negligible for n = 0, . . . , 3 when r > S . We conclude that the non-linear convective term can be neglected in the outer problem, because it is sub-leading within and beyond the matching region for n = 0, . . . , 3. In the inner problem, in general, the non-linear convective term cannot be neglected for n ≤ 3.

C. Matching

The outer solution [START_REF] Rosén | Numerical analysis of the angular motion of a neutrally buoyant spheroid in shear flow at small Reynolds numbers[END_REF] provides boundary conditions for the inner problem in the matching region, at r ∼ 1/ . The leading-order terms of u in,i at large r must match Eq. [START_REF] Rosén | Numerical analysis of the angular motion of a neutrally buoyant spheroid in shear flow at small Reynolds numbers[END_REF]. We now discuss the matching order by order. To leading order 0 the inner problem reads:

-∂ i p (0) in +∂ j ∂ j u (0) in,i = 0 . ( 18
)
The boundary and matching conditions are given by:

u (0) in,i = ε ijk ω (0) j r k -u (∞) i for r ∈ S , (19a)
u (0) in,i ∼ T i (0) as r → ∞ . ( 19b
)
From Eq. (19b) we see that T i (0) is the leading-order contribution to the Stokes problem of a spheroid freely rotating in a simple shear with its symmetry axis aligned with the ê3 -axis. The matching condition (19b) determines the dipole strength D introduced in the previous Section. One finds D as a function of the aspect ratio λ:

D = 16π(λ 2 -1) 3 3λ 3 [5λ -7λ 3 + 2λ 5 + 3 √ λ 2 -1 arccosh(λ)] (20a)
for a prolate spheroid (λ > 1), and

D = - 16π(1 -λ 2 ) 3 3[5λ -7λ 3 + 2λ 5 -3 √ 1 -λ 2 arccos(λ)] (20b)
for an oblate spheroid (λ < 1). The corresponding angular velocity is the Jeffery-result ω (0) = -1 2 [START_REF] Jeffery | The motion of ellipsoidal particles immersed in a viscous fluid[END_REF]. The order-1 problem reads:

-∂ i p (1) in +∂ j ∂ j u (1) in,i = 0 , (21a) u (1)
in,i = ε ijk ω

(1)

j r k for r ∈ S , (21b)
together with the matching condition

u (1) in,i ∼ T (1) i as r → ∞ . ( 21c
)
Since T

(1) i vanishes [see Eq. (15b)] there is no term in Eq. ( 14) to which the flow u [START_REF] Einarsson | Rotation of a spheroid in a simple shear at small Reynolds number[END_REF] in,i produced by the rotating particle can be matched. We conclude that u [START_REF] Einarsson | Rotation of a spheroid in a simple shear at small Reynolds number[END_REF] in,i ≡ 0. It follows that the 1 -contribution to the angular velocity must vanish, ω (1) = 0.

The second-order problem is inhomogeneous:

-∂ i p (2) in +∂ j ∂ j u (2) in,i = A (∞) jk r k ∂ j u (0) in,i +A (∞) ij u (0) in,i +u (0) in,j ∂ j u (0) in,i , (22a)
with boundary and matching conditions:

u (2) in,i = ε ijk ω (2) j r k for r ∈ S , (22b) u
(2) in,i ∼ T i (2) as r → ∞ . ( 22c
)
The solution is the sum of a homogeneous and a particular part:

u (2) in,i = u (2) h,i + u (2)
p,i . The particular solution is found to be u 2) +O(1/r 2 ). The homogeneous part solves the Stokes problem

(2) p,i = T i (
-∂ i p (2) h +∂ j ∂ j u (2) h,i = 0 (23a)
with boundary conditions

u (2) h,i = ε ijk ω (2) j r k -u (2)
p,i for r ∈ S .

The explicit solution of the homogeneous problem shows that u

h,i is asymptotic to O(1/r 2 ) as r → ∞. This implies that the matching condition (22c) is fulfilled by u (2) p,i alone.

Direct computation confirms that the 2 -contribution to the angular velocity of the particle must vanish for torque-free rotation, ω (2) = 0. We remark that there are, however, examples where there is a linear inertial correction to the angular velocity [START_REF] Candelier | Angular dynamics of small particles in turbulence[END_REF].

So far we have found that ω (0) = -1 2 , and ω (1) = ω (2) = 0. Now consider the third order in . At this order the problem is homogeneous since u [START_REF] Einarsson | Rotation of a spheroid in a simple shear at small Reynolds number[END_REF] in,i ≡ 0. The equations to solve are

-∂ i p (3) in +∂ j ∂ j u (3) in,i = 0 , (24a) u (3)
in,i = ε ijk ω

(3)

j r k for r ∈ S , (24b)
together with the matching condition

u (3) in,i ∼ T i (3) as r → ∞ . (24c)
Eq. (15d) shows that T

(3) i is a linear flow. The thirdorder problem (24) is thus equivalent to that of a freely rotating speroid in the linear flow A ij r j as → 0. The solution of this problem is known, it is discussed in the next Section.

D. Angular velocity

The third-order correction to the angular velocity is computed by requiring that the third-order torque τ (3) ≡ τ (3) ê3 vanishes. The third-order inner problem is equivalent to the problem of determining the torque on a freely rotating spheroid in a linear flow in the creeping-flow limit. The solution is known (see e.g. p.64 in [START_REF] Kim | Microhydrodynamics: principles and selected applications[END_REF]). In the log-rolling orbit we have:

τ (3) = 8π 3D 20π (A 21 -A 12 ) -ω (3) . ( 25
)
Setting Eq. ( 25) to zero yields ω (3) , and thus:

ω ∼ -1 2 + 3D 20π (A 21 -A 12 ) 3 (26)
in dimensionless variables. In Section IV we show how to evaluate the constant coefficients A ij . The result is

A 12 = 0.0328 , A 21 = 0.1408 . ( 27
)
It follows that

ω∼ -1 2 + 0.0540 3D 10π 3 . (28)
In dimensional variables this result corresponds to Eq. ( 1).

IV. SOLUTION OF THE OUTER PROBLEM

In this Section we determine the solution of the outer problem (13) in the limit → 0. The results are Eqs. [START_REF] Rosén | Numerical analysis of the angular motion of a neutrally buoyant spheroid in shear flow at small Reynolds numbers[END_REF].

Since Eq. ( 13) is linear it can be solved by Fourier transform (we employ the symmetric convention). Using incompressibility to eliminate the pressure, we find for the Fourier transform û o,i (k):

2 k 1 ∂ û o,i ∂k 2 -δ 1i - 2k 1 k i k 2 û o,2 (29)
= k 2 û o,i + iD (2π) 3/2 k 2 2k i k 1 k 2 -k 2 (k 1 δ i2 + k 2 δ i1 ) .
To determine the terms T (n) i in Eq. ( 14) we expand û o,i in . We see in the next Section, however, that a regular expansion does not give all required terms. There are terms involving generalised functions that must be determined separately (Section IV B).

A. Even-order terms

It follows from Eq. ( 29) that the leading order in the -expansion of û o,2 takes the form:

û o,2 = - iD (2π) 3/2 k 1 k 4 2k 2 2 -k 2 + O( ) . ( 30
)
This determines T (0) 2 (k):

T (0) 2 (k) = iD (2π) 3/2 k 1 k 4 k 2 -2k 2 2 . ( 31a
)
This is the Fourier transform of Eq. (15a). The next term in the expansion reads:

T (2) 2 (k) = - 4 iD (2π) 3/2 k 2 1 k 2 k 8 k 2 -k 2 2 . ( 31b
)
This is the Fourier transform of Eq. (15c). The next order in the expansion is:

T (4) 2 (k) = - 4 iD (2π) 3/2 k 3 1 k 12 k 2 -6k 2 2 . ( 31c
)
In a similar way we expand u o,1 . This yields

T (0) 1 (k) = iD (2π) 3/2 k 2 k 4 (k 2 -2k 2 1 ) , (32a)
T (2) 1 (k) = 4 iD (2π) 3/2 k 3 1 k 2 2 k 8 , (32b)
T (4) 1 (k) = 4 iD (2π) 3/2 k 2 1 k 2 k 12 k 4 1 +(2k 2 3 -5k 2 2 )k 2 1 +k 2 3 (k 2 2 +k 2 3 ) . (32c)
We see that a regular expansion of Eq. ( 30) yields only even orders T (2n) i

. Yet odd orders T (2n-1) i are required to match the inner solution. Terms corresponding to odd powers in are computed in the following Sections.

B. Odd-order terms

Odd terms in the expansion [START_REF] Subramanian | The influence of the inertially dominated outer region on the rheology of a dilute dispersion of low-reynoldsnumber drops or rigid particles[END_REF] are not captured by the expansion in described in the previous Section because the odd-order terms are zero everywhere except at isolated points in k-space. To compute the odd-order terms we subtract the lowest term(s) in the expansion [START_REF] Subramanian | The influence of the inertially dominated outer region on the rheology of a dilute dispersion of low-reynoldsnumber drops or rigid particles[END_REF] to make the next odd order leading. This new leading order is then found by taking the limit → 0.

This idea was pioneered by Childress [START_REF] Childress | The slow motion of a sphere in a rotating, viscous fluid[END_REF] and Saffman [START_REF] Saffman | The lift on a small sphere in a slow shear flow[END_REF], and employed by Lin et al. [START_REF] Lin | Simple shear flow around a rigid sphere: inertial effects and suspension rheology[END_REF]. These papers start from an ansatz of the outer solution written in terms of a 'stretched' configuration-space variable r ≡ r, see for example Eq. (3.12) in Ref. [START_REF] Saffman | The lift on a small sphere in a slow shear flow[END_REF]. Here we pursue an alternative path that does not absorb the into a new variable r and proceeds with generalised functions in k-space [START_REF] Candelier | Note on the method of matched-asymptotic expansions for determining the force acting on a particle[END_REF]. We define:

û(q) o,i (k) ≡ 1 q-1 û o,i (k) - q-1 k=0 k T (k) i (k) (33) ∼ T (q) i (k) + 2 T (q+1) i (k) + . . . .
To extract an odd term corresponding to q = 2n -1 one can take the limit → 0 [6]:

T (2n-1) i (k) = lim →0 -1 û(2n-1) o,i (k) . ( 34a
)
It turns out, however, that an alternative procedure has distinct advantages:

T (2n-1) i (k) = lim →0 d d -1 û(2n-2) o,i (k) . (34b)
Eqs. (34) are evaluated as follows. We derive a differential equation for û(

by substituting Eq. (33) into Eq. (29):

2 k 1 ∂ û(2n-1) o,i ∂k 2 -k 2 û(2n-1) o,i -2 δ 1i - 2k 1 k i k 2 û(2n-1) o,2 (35)
= -k 2 T (2n) i .

For i = 1, 2 the solutions of Eq. ( 35) are:

û(2n-1) o,2 (k) = - 1 k 2 ∞ 0 dξ e -ξ k 2 1 ξ 2 /3+k1k2ξ+k 2 / 2 k 4 T (2n) 2 (k) k2=k1ξ+k2 , (36a) û(2n-1) o,1 (k)
= - ∞ 0 dζ e -ζ k 2 1 ζ 2 /3+k1k2ζ+k 2 / 2 k 2 T (2n) 1 (k) + k 2 -2k 2 1 k 1 k 2 û(2n-1) o,2 (k) k2=k1ζ+k2 , (36b)
Then, to determine T (2n-1) i we must take the limit → 0 in Eqs. (34). To this end we use the fact that û(2n-1)

o,i (k
) are integrable functions of k. For such a function the limit → 0 can be expressed in terms of generalised functions. We need the following three relations. First,

lim →0 1 3 f (k/ ) = A δ(k) with A = R 3 d 3 k f (k) . (37a)
Second, consider a function g that integrates to zero,

R 3 d 3 k g(k) = 0 . (37b)
In this case we use

lim →0 g(k/ ) 4 = A i ∂δ(k) ∂k i with A i = - R 3 d 3 k k i g(k) . (37c)
Third, we make use of the identity:

lim →0 d d h(k/ ) 3 = A i ∂δ(k) ∂k i with A i = - R 3 d 3 k k i h(k) .
(37d) This relation follows from Eq. (37a) for a differentiable function h(k) (Appendix A). To make use of Eqs. (37) we first derive a scaling relation for T

i applying the results of Section IV A:

T (2n) i (k) = 1 2n+1 T (2n) i (k/ ) . (38)
This relation together with (36) implies:

û(2n-1) o,i (k) = 1 2n-1 û(2n-1) o,i (k/ ) . (39)
Eq. (39) allows us to bring (34) into a form where Eqs. (37) can be directly applied. This yields:

T (2n-1) i (k) = lim →0 1 2n û(2n-1) o,i (k/ ) , (40a)
T (2n-1) i (k) = lim →0 d d 1 2n-1 û(2n-2) o,i (k/ ) , (40b)
Substituting n = 1 into (40a) gives:

T (1) i (k) = lim →0 1 2 û(1) o,i (k/ ) . (41)
This limit evaluates to zero, and we conclude that

T (1) i (r) ≡ 0 [Eq. (15b)]. Now consider T (3) i (k).
First, we use the fact that T (1) i (k) = 0 together with Eq. ( 33) to find the relation

û(2) o,i (k) = -1 û(1) o,i (k) . (42)
It follows from Eqs. ( 40) and ( 42) that:

T (3) i (k) = lim →0 1 4 û(3) o,i (k/ ) , (43a)
T (3) i (k) = lim →0 d d 1 3 û(1) o,i (k/ ) , (43b)
Finally, applying Eqs. (37c) and (37d) we find

T (3) i (k) = A ij ∂ ∂k j δ(k) , (44)
together with the two equivalent expressions for the normalisation

A ij A ij = - R 3 d 3 k k j û(3) o,i (k) , (45a)
A ij = - R 3 d 3 k k j û(1) o,i (k) . (45b)
The Fourier transform of Eq. ( 44) is Eq. (15d), with

A ij = -i (2π) 3/2 A ij . (46)
The prefactor comes from the Fourier transform (we use the symmetric convention). The integrals A 21 and A 12 must be evaluated numerically. It is interesting to note that Eq. (45a) agrees with Eq. (3.37) in Ref. [START_REF] Lin | Simple shear flow around a rigid sphere: inertial effects and suspension rheology[END_REF] (Lin et al. use a different convention for the Fourier transform). This confirms that our expansion in k-space yields the same result for a sphere as that obtained by Lin et al., using Saffman's asymptotic matching method. We could not perform the numerical integrations of Eq. (45a) with sufficient accuracy (we estimate the relative error of our numerical integrations to be of the order of ten percent). We therefore suspect that the error in the coefficient computed by Lin et al. occurred in the numerical integration. On the other hand, substituting Eq. (45b) into Eq. ( 26) and we find an expression that, for a spherical particle, is equivalent to the result obtained by by Stone, Lovalenti & Brady [START_REF] Stone | Inertial effects on the rheology of suspensions and on the motion of individual particles[END_REF]. We emphasise that Eq. (45b) expresses the coefficients A ij (and thus the angular velocity) in terms of û(1) o,i [Eq. ( 33)]. It is not necessary to compute û(3) o,i , as suggested by Saffman [START_REF] Saffman | The lift on a small sphere in a slow shear flow[END_REF] and Lin et al. [START_REF] Lin | Simple shear flow around a rigid sphere: inertial effects and suspension rheology[END_REF], and carried out in Eq. (45a). This detail makes a distinct difference. The integrals in Eq. (45b) are much simpler to evaluate numerically. We find:

A 12 = i 0.517 and A 21 = i 2.220 , (47)
Using Eq. ( 46) we obtain the coefficients quoted in Section III D.

V. DIRECT NUMERICAL SIMULATIONS

The numerical simulations are performed using the commercial finite-element software package Comsol Multiphysics 5.1. The simulation domain is a cube with side length L. Since κ = 2a/L, the side length is 2κ -1 in dimensionless variables. We impose velocity boundary conditions in the shear direction, u 1 = ±κ -1 at r 2 = ±κ -1 , and periodic boundary conditions in the flow and vorticity directions at r 1,3 = ±κ -1 . On the particle surface no-slip boundary conditions are used. For given values of Re s and κ we obtain the steady-state solution of Eq. ( 3), together with the corresponding angular velocity ω.

The numerical mesh consists of tetrahedral elements. The spatial resolution changes, from ∆r min close to the particle to ∆r max far from the particle. All simulations are performed using ∆r min = 0.06 and ∆r max = 0.1κ -1 . To assess the numerical precision we also perform simulations with ∆r min = 0.12. The difference |ω(∆r min = 0.12) -ω(∆r min = 0.06)| estimates the absolute numerical error of the angular velocity. We estimate this error to be smaller than 10 -6 in dimensionless units.

We compute the creeping-flow limit ω (0) of the angular velocity in the finite system by evaluating the angular velocity for a fluid with zero mass density. We estimate that the absolute error of the finite-size corrections is smaller than 10 -6 in dimensionless units.

VI. DISCUSSION AND CONCLUSIONS

We computed the inertial corrections to the angular velocity of a spheroid log-rolling in an unbounded simple shear to order 3 = Re 3/2 s . Our main result, Eq. ( 1), is shown as a function of the particle aspect ratio λ by the solid red line in Fig. 3. We observe excellent agreement with the results of our direct numerical simulations (shown as open symbols). Furthermore, we see that the inertial correction is largest for the spherical particle. This reflects the fact that the surface area |S | of the particle is largest when λ = 1. The surface area is also shown in Fig. 3, as a black dashed line. As λ → ∞ the inertial correction vanishes. This is expected because the latter limit corresponds to slender rods, which do not couple back to the fluid through the boundary conditions. For extended bodies, fluid elements near the rotating particle surface are centrifuged away from the rotation axis due to fluid inertia, thereby slowing the particle down. In the slender-disk limit (λ → 0) the inertial correction approaches a finite constant ∼ Re 3/2 s . The results of the numerical simulations for a spherical particle are plotted as a function of Re s in Fig. 2, for three different system sizes: κ = 2a/L = 0.01, 0.025, and 0.05. Shown is the difference between ω and the value of the angular velocity in creeping-flow limit, ω (0) . The inset of Fig. 2 shows ω (0) as a function of κ. A fit to

ω (0) (κ) = -1 2 + Cκ 3 , (48)
gives C ≈ 0.22. For Re s between 5×10 -4 and 5×10 -2 the numerical data for the largest system agree well with the theoretical prediction (28). This conclusion is supported also by the numerical results reported in Ref. [START_REF] Mikulencak | Stationary shear flow around fixed and free bodies at finite Reynolds number[END_REF].

But at smaller values of the shear Reynolds number the numerical data lies below the theory. The reason is that L < S at very small values of Re s . In this case the outer region does not exist and the problem becomes a regular perturbation problem. We therefore expect that the inertial correction is quadratic in this regime:

ω = ω (0) (κ) + b(κ)Re 2 s when Re s κ 2 . ( 49
)
The thin solid (blue) lines in Fig. 2 show that this is indeed the case. Fitting Eq. (49) to the numerical data we find that the coefficient b(κ) increases monotonously as κ decreases. The numerical values of b(κ) are given in Table I. Also shown in Fig. 2 are experimental results of Poe & Acrivos [START_REF] Poe | Closed-streamline flows past rotating single cylinders and spheres: inertia effects[END_REF] for a sphere in a shear flow at shear Reynolds numbers of order unity. At these Reynolds numbers the theory (1) for λ = 1 is no longer valid. The results of our numerical simulations are slightly below the experimental results. We do not know the reason for this discrepancy, but note that Poe & Acrivos [START_REF] Poe | Closed-streamline flows past rotating single cylinders and spheres: inertia effects[END_REF] comment on the fact that a re-circulation flow observed in their experiments at Re s ≈ 5 may have caused a spurious increase of the angular velocity.

In summary, we observe good agreement between the numerical simulations and theory where expected, and we find reasonable agreement between the simulations and the experiment at Re s of order unity.

The details of our matching method differ from the standard method [START_REF] Proudman | Expansions at small Reynolds numbers for the flow past a sphere and circular cylinder[END_REF][START_REF] Saffman | The lift on a small sphere in a slow shear flow[END_REF]. We conclude by briefly summarising the advantages of our approach. We have described two different schemes of calculating the thirdorder contribution to the angular velocity. The first scheme (Section IV B) expresses the angular velocity in terms of û(3) o,i (k) [Eq. ( 33)]. This approach is equivalent to the standard asymptotic matching scheme used by Lin et al. [START_REF] Lin | Simple shear flow around a rigid sphere: inertial effects and suspension rheology[END_REF]. It results in integral expressions that we could not evaluate accurately enough.

The second scheme allows us to express the third-order contribution to the angular velocity in terms of û(1) o,i (k). This leads to integrals that are much easier to evaluate. Saffman [START_REF] Saffman | The lift on a small sphere in a slow shear flow[END_REF] remarks that the calculation of such terms in this outer expansion is "a matter of great difficulty". Our new scheme is therefore of practical importance when performing asymptotic matching, as we avoid computing two orders in = √ Re s altogether. Our approach shares this advantage with an alternative method that is based on the reciprocal theorem [START_REF] Stone | Inertial effects on the rheology of suspensions and on the motion of individual particles[END_REF], and leads to identical integral expressions. This provides provides an interesting link between the standard asymptotic matching method and approaches based on the reciprocal theorem. But it remains to be investigated how general this correspondence is. We expect that our scheme offers advantages in different problems too. One example is the evaluation of the viscosity and the normal stresses in a dilute suspension of spheroids. For spherical particles these were calculated in Refs. [START_REF] Lin | Simple shear flow around a rigid sphere: inertial effects and suspension rheology[END_REF][START_REF] Stone | Inertial effects on the rheology of suspensions and on the motion of individual particles[END_REF][START_REF] Subramanian | The influence of the inertially dominated outer region on the rheology of a dilute dispersion of low-reynoldsnumber drops or rigid particles[END_REF].

In an unbounded shear the log-rolling orbit is stable for oblate particles, but unstable for prolate particles. Direct numerical simulations [START_REF] Rosén | Numerical analysis of the angular motion of a neutrally buoyant spheroid in shear flow at small Reynolds numbers[END_REF] show that the corresponding stability exponent of the log-rolling orbit has a Re 3/2 s -correction which is of the same order as the coefficients of the linear contribution. It would therefore be of interest to compute this correction from first principles. It is also important to mention that the direct numerical simulations are performed for a finite system, a bounded shear. It would be of great interest to generalise the results obtained here, but also those reported in Refs. [START_REF] Einarsson | Rotation of a spheroid in a simple shear at small Reynolds number[END_REF][START_REF] Einarsson | Effect of weak fluid inertia upon Jeffery orbits[END_REF][START_REF] Candelier | The role of inertia for the rotation of a nearly spherical particle in a general linear flow[END_REF], to finite systems. The first step is to compute the small-Re 2 s asymptotes in Fig. 2.

  FIG. 1.Spheroid log rolling in a simple shear. The r2coordinate is aligned with the shear direction, and the r3coordinate is the negative vorticity direction. The axis of symmetry of the particle and its angular velocity ω = ωê3 are aligned with the flow vorticity.

  FIG. 3. Re3/2

TABLE I .

 I Numerical values of the parameter b(κ) in Eq. (49), obtained by fitting the results of the direct numerical simulations. For κ = 0.01 the fit used values of Res smaller than 10 -3 . For the other values of κ the data was fitted up to Res = 2.5 × 10 -3 .

	κ	0.01 0.025 0.03 0.04 0.05 0.06 0.08 0.1
	b(κ) 1.48 0.60 0.49 0.35 0.27 0.21 0.14 0.10

(∂ t u i +u j ∂ j u i ) = -∂ i p + ∂ j ∂ j u i , ∂ j u j = 0 ,(3)

ACKNOWLEDGMENTS

We thank Howard Stone for discussions, and for making the preprint [5] available to us. We thank Kirsten Mehlig for help with the analysis of the numerical data. JM, JE and BM acknowledge support by Vetenskapsrådet [grant number 2013-3992], Formas [grant number 2014-585], by the grant 'Bottlenecks for particle growth in turbulent aerosols' from the Knut and Alice Wallenberg Foundation, Dnr. KAW 2014.0048, and by the MPNS COST Action MP1305 'Flowing matter'. TR and FL acknowledge financial support from the Wallenberg Wood Science Center, Bengt Ingeström's Foundation, and ÅForsk. The computer simulations were performed using resources provided by the Swedish National Infrastructure for Computing at the National Supercomputer Center in Sweden.

Appendix A

In this appendix we derive relation (37d). To this end we consider the action of lim →0

This expression can written in the form of a divergence:

We use Eq. (37a) to evaluate the last expression further:

In summary we find:

with

This is Eq. (37d).