# Dynamics of spheroids in pressure driven flows of shear thinning fluids

Vishal Anand and Vivek Narsimhan 111To whom the correspondence should be
addressed. mailto:vnarsim<EMAIL_ADDRESS>
Davidson School of Chemical Engineering, Purdue University,
West Lafayette, Indiana 47907, USA

###### Abstract

Particles in inertialess flows of shear thinning fluids are a model
representation for several systems in biology, ecology, and microfluidics. In
this paper, we analyse the motion of a spheroid in a pressure driven flow of a
shear thinning fluid. The shear thinning rheology is characterised by the
Carreau model. We use a combination of perturbative techniques and the
reciprocal theorem to delineate the kinematics of prolate and oblate
spheroids. There are two perturbative strategies adopted, one near the zero
shear Newtonian plateau and the other near the infinite shear Newtonian
plateau. In both limits, we find that a reduction in effective viscosity
decreases the spheroid’s rotational time period in pressure driven flows. The
extent to which shear thinning alters the kinematics is a function of the
particle shape. For a prolate particle, the effect of shear thinning is most
prominent when the spheroid projector is aligned in the direction of the
velocity gradient, while for an oblate particle the effect is most prominent
when the projector is aligned along the flow direction. Lastly, we compare the
tumbling behavior of spheroids in pressure driven flow to those in simple
shear flow. While the time period decreases monotonically with Carreau number
for pressure driven flows, the trend is non monotonic for shear flows where
time period first increases at low Carreau number and then decreases at high
Carreau numbers. Shear thinning does not resolve the degeneracy of Jefferey’s
orbits.

##  1 Introduction

Rigid, orientable particles experience a bevy of interesting phenomena when
placed in highly viscous, structureless fluids [1]. For example, a single rod
sediments at the same initial orientation in such fluids, while a single rod
tumbles in a periodic orbit (Jeffrey orbit) in shear flow. The reason for
these behaviors arise from the symmetry, linearity, and reversibility of the
Stokes equations [2, 3, 4]. When the assumptions of the Stokes equations break
down – i.e., when fluid inertia is present or the fluid exhibits non-Newtonian
rheology –the above behavior may no longer be valid  [5] .

By now, there are plenty of studies that investigate how different particle
shapes (e.g., fibers, sheets, spheres, spheroids) move under sedimentation,
linear flows, and quadratic flows when the fluid has small but non-zero
inertia. For example, a fibre-like particle was analysed in shear flow of a
Newtonian fluid with weak inertia using a slender body approximation and the
reciprocal theorem [6]. It was shown that the fibre drifts towards the
velocity-shear gradient plane due to inertia and then stops rotating
altogether above a critical Reynolds number. A similar analysis was performed
for spheroids in simple shear flow [7]. The analysis was later extended to
general linear flows, which are a combination of extensional and rotational
flows [8]. The analyses showed that both fluid and particle inertia forced the
prolate particle to tumble in velocity -velocity gradient plane, while the
oblate particle was forced into a log rolling motion. A linear stability
analysis validated the stability of these orientations for prolate and oblate
particles [9, 10].

When the suspending fluid is no longer Newtonian, the fluid can experience
phenomena such as normal stress differences, shear thinning, extensional
thickening, and/or time-dependent viscoelasticity. Each of these effects play
a major role in altering the motion of rigid, orientable particles. We will
not provide a complete survey on the microhydrodynamics of such particles in
non-Newtonian fluids – instead, we refer the readers to treatises [11, 12] and
provide a brief summary here. Generally, much of the effort on the past few
decades has focused on the effect of normal stresses on particle motion [13,
14, 15]. Normal stresses give rise to two different phenomena in low Reynolds
number flows – namely cross stream migration and steady state orientation.
Cross-stream migration refers to the lift forces on particles flowing in
channel. This phenomenon occurs when there is a gradient in normal stresses,
and thus is generally observed when there is a gradient in strain rates –
i.e., a quadratic flow [16]. During sedimentation and/or shear flow, normal
stresses can also give rise to a hydrodynamic torque on orientable particles.
This situation leads to a stable orientation [16, 17]. The analytical solution
for the sedimentation of spheroids in a quiescent background of a second order
fluid was furnished by [18], while those for linear flows by [19]. Recently,
these studies were extended by [20], who furnished analytical formulae for the
polymeric force and torque on a spheroid, in quadratic flows of second order
fluids in the limit of small Weissenberg number. A recent survey of the
research into viscoelastic microhydrodynamics is furnished in [12].

There has been interest in the microhydrodynamics of particles in fluids with
variable viscosity, where the variation in viscosity is either spatial [21,
22, 23] or due to the shear thinning rheology of the fluid [24, 25]. Out of
the two, the focus in this paper is the motion of spheroids in shear thinning
fluids; such systems are prevalent in biology, microfluidics, and ecology [26,
27]. There particularly has been work done for microswimming (active
particles)[28], due to the intuitive but naive speculation that a reduction in
viscosity afforded by shear thinning allows the swimmer to swim faster without
expending as much as energy [29, 27]. For spherical squirmers, an asymptotic
analysis revealed both enhancement and reduction of swimming velocity
depending on the surface actuation of swimmers [30, 31]. There are shear rates
at which the swimming velocity gets optimized [31]. Thus, the naive hypothesis
that shear thinning enhances the swimming speed is faulty and must be used
with caution [31, 30]. Very recently, the swimming characteristics of
Purcell’s swimmer in a shear thinning fluid was analysed where it was revealed
that unequal arm rotation rates induce a net vertical displacement, which is
not present in Newtonian fluids [32].

Despite the recent progress on particle dynamics in shear thinning fluids, it
is clear that this area of research is relatively underexplored with many
challenging problems remaining to be investigated. For instance, most of the
asymptotic analyses of particle dynamics explore the Carreau model of fluid
rheology, performing a perturbation expansion in the small Carreau number
(i.e., low shear rate) limit. We still do not know what happens when the
Carreau number is not small. Moreover, for passive particles, the analysis has
been restricted to spheres, with the notable exception of [25]. Other particle
shapes, like oblate spheroids, have not been investigated with respect to
their interaction with shear thinning flows. Finally, we understand that
microhydrodynamics of particles during sedimentation [24] and linear shear
flows [25] of shear thinning fluids has been analysed, but more complicated
flows like pressure driven flows have not been discussed. In this context, we
propose to answer the following research questions in this paper: How does the
orientational kinematics of prolate and oblate spheroids vary in pressure
driven flows of shear thinning fluids in both the small and large Carreau
number limits? How is this behavior different than in the simple shear flow
case?

The paper is organized as follows. First we introduce the problem formally in
Sec. 2. The definition of prolate and oblate spheroids are discussed here as
well as the fluid rheological model (Carreau model). Next, in Sec. 3, we show
how the pressure driven flow is altered due to shear thinning alone, without
the presence of any particle. This section also introduces the particle into
this flow field, and uses the reciprocal theorem to quantify the shear
thinning correction to the particle kinematics. In the next section, we
discuss the algorithm of the code (Sec. 4), which we use to obtain the results
in Sec. 5. In Sec. 6, we present some results pertaining to the tumbling time
period of spheroids in simple shear flows in both the small and large Carreau
number regimes and discuss the difference between the trends observed in shear
flows vis-a-vis those observed in pressure driven flows. A discussion and
conclusion follow in Sec. 7.

##  2 Problem Statement

###  2.1 Problem geometry and particle definition

The system under investigation is a spheroid in a pressure driven flow of a
shear thinning fluid. The flow is steady and inertialess, and the particle is
neutrally buoyant and passive. Figure 1a illustrates the geometry of the
system. We investigate a slit-like channel of length $l$, height $2h$, and
infinite width. The coordinate system is positioned such that the $x$-axis
aligns with the midplane of the channel while the $y$-axis aligns with the
height of the channel. A pressure $p=\Delta P$ is imposed at the inlet of the
channel $x=0$, while the outlet pressure at $x=l$ is $p=0$.

The spheroid will start at a position $(x_{0},y_{0},z_{0})$ in the channel.
The lengths of the three semi-axes are $a,b,$ and $c$, where $b=c$. The
unequal axis – i.e., the direction of the $a$-axis –- is known as the
projector. A prolate spheroid has its unequal axis the longest ($a>b$), while
an oblate spheroid has its unequal axis the shortest ($a<b$). Apart from the
semi-axes ($a,b,c$), it is also possible to parameterize the spheroid shape
using two other quantities $R$ and $A_{R}$. $R$ is the radius of the
equivalent sphere with the same volume as the spheroid, which means
$V=\frac{4}{3}\pi R^{3}=\frac{4}{3}\pi abc$. The quantity $A_{R}$ is the ratio
of the projector axis to the other two axes: $A_{R}=a/b$. By this definition,
prolate spheroids have $A_{R}>1$, while oblate spheroids have $A_{R}<1$. The
two systems of particle parameterization are connected by

$a=RA_{R}^{2/3},\quad b=c=RA_{R}^{-1/3}$ (1)

The orientation of the spheroid is characterized by the ordered pair
($\theta,\phi$). Here $\theta\in(0,\pi)$ is the angle from the $z$-axis, also
denoted as the co-latitude angle or the polar angle, while $\phi\in(0,2\pi)$
is the angle in $x-y$ plane from the $x$-axis, also known as the azimuth
angle. Figure 2 summarizes the geometrical definitions for the spheroids.

We will solve for the motion of the spheroid in the limit when the particle
size is much smaller than the channel height – i.e, $R\ll h$ . In this
situation, one can neglect the hydrodynamic interaction with the wall and
treat the particle as if it were in an unbound fluid with a background
velocity $u_{i}^{\infty}$ given by the flow field in the channel. This
approximation incurs an error of $O((R/h)^{3})$ for the rigid body motion (see
Kim and Karilla [4, Ch. 12]). Thus, in this problem, we will often switch
between two coordinate systems. When solving the translational and rotational
velocity of the particle, we will use Figure 1b and treat the particle in an
unbound medium with the origin at the particle’s center of mass and local
coordinate system aligned with the particle’s semi-axes. Once we determine
these quantities, we will go back to the channel coordinates (Figure 1a) and
update the particle position and orientation over a time step $\Delta t$. The
process will repeat until we track the particle motion for a given time
interval. In this context, the microhydrodynamic problem encountered here may
be considered as a case of fluid structure interaction problem [33, 34, 35,
36, 37, 38]; the flow field exerts pressure/force on the structure leading to
the motion (or deformation, if the solid is soft) of the particle. The next
subsection describes the rheological model for the fluid in the channel.

(a)

(b)

Figure 1: System geometry in a shear thinning fluid. (a) Lab coordinate
system. The channel has half width $h$ and length $l$, and the coordinate axes
are placed at the center of the channel inlet. A spheroid is initially placed
at location $(x_{0},y_{0},z_{0})$ and its position and orientation is tracked
over time. The dimensions are not to scale. (b) Particle coordinate system.
When solving for the particle’s rigid body motion, we assume the particle is
in an unbound fluid with a velocity $u_{i}^{\infty}$ given by the flow in the
channel. The coordinates align with the particle’s semi-axes and the origin is
the particle’s center of mass.

(a)

(b)

Figure 2: Definition of the particle shape and orientation for both (a)
prolate and (b) oblate spheroids. The projector makes an angle $\phi$
(azimuth) in the $x-y$ plane with flow direction $x$, and an angle $\theta$
(polar) with the $z$ direction. For both prolate and oblate spheroids, the
relationship is also shown between the two different systems of
parameterizations, based on $a,b,c$ and on $R$ and $A_{R}$ respectively.

###  2.2 Fluid rheology

####  2.2.1 Carreau-Yasuda model

The background fluid is complex [39, 40, 41]. Specifically, the fluid is
shear-thinning – i.e., the fluid exhibits a reduction in the apparent
viscosity with increase in the applied shear stress. Such fluids are also
called pseudoplastics and common examples include polymers, blood, and ketchup
among others [42](see also [43, 44]). The constitutive equation for the stress
is the same as a Newtonian fluid, except that the viscosity is a function of
the strain rate:

$\tau_{ij}=\eta(\dot{\gamma})\dot{\gamma}_{ij}–p\delta_{ij};$ (2)

In the above equation, $p$ is the pressure, $\dot{\gamma}_{ij}=\frac{\partial
u_{i}}{\partial x_{j}}+\frac{\partial u_{j}}{\partial x_{i}}$ is twice the
rate of strain tensor, and $\eta$ is the viscosity. The viscosity is a
function of the magnitude of strain rate, given by
$\dot{\gamma}=\sqrt{\frac{1}{2}\dot{\gamma}_{ij}\dot{\gamma}_{ij}}$, where
Einstein convention is assumed.

Rheological models for shear thinning viscosity include the power law model of
Ostwald–de Waele, the Herschel-Bulkley model, and the Carreau-Yasuda model
[39]. Here, we employ the Carreau-Yasuda model, mainly because this model
captures both Newtonian plateaus at the beginning and the end of shear
thinning regimes. The constitutive equation for the Carreau-Yasuda model is
given as:

$\eta=\eta_{\infty}+(\eta_{0}-\eta_{\infty})[1+\lambda_{t}^{2}\dot{\gamma}^{2}]^{\frac{n-1}{2}}$
(3)

In this model, the fluid behaves Newtonian at low and high shear rates, with a
power-law region in-between. The quantities $\eta_{0}$ and $\eta_{\infty}$ are
the zero shear rate and infinite shear rate viscosities, respectively. The
index $n<1$ determines the rate of decay in the power-law region, while
$\lambda_{t}$ is the time constant that determines the strain rate at which
the power-law region occurs.

####  2.2.2 Non-dimensional form of Carreau-Yasuda model

The dimensionless form of the Carreau-Yasuda model is given by:

$\eta=1-\epsilon+\epsilon[1+Cu^{2}\dot{\gamma}^{2}]^{\frac{n-1}{2}}$ (4)

where the viscosity is rendered dimensionless by $\eta_{0}$ and the shear rate
$\dot{\gamma}$ is rendered dimensionless by a characteristic shear rate
$\dot{\gamma}_{c}$. The characteristic shear rate is determined by the problem
under consideration (such as shear flow or pressure driven flow). There are
two new dimensionless numbers introduced above. The relative viscosity drop

$\epsilon=1-\beta=\frac{\eta_{0}-\eta_{\infty}}{\eta_{0}}$ (5)

represents the fractional drop in viscosity between the zero and infinite
shear rate limits. The Carreau number

$Cu=\lambda_{t}\dot{\gamma_{c}}$ (6)

describes the shear-thinning regime at the characteristic strain rate
$\dot{\gamma}_{c}$. When $Cu\ll 1$ or $Cu\gg 1$, the characteristic shear rate
is located near the Newtonian plateaus, while for $Cu\sim O(1)$ the
characteristic shear rate is in the power-law regime. A summary of the
dimensionless numbers in the Carreau model is shown in Table 1. For
illustration, typical values of the rheological parameters of Carreau model
pertaining to aqueous solution of $0.3\%$ xanthan gum (molar mass = $933.748$
g/mol) are $n=0.402$ , $\lambda_{t}=239$ s and $\beta=1.35\times 10^{-4}$ [45]
(see also [46]).

Table 1: Summary of dimensionless rheological parameters for Carreau shear thinning fluid. Dimensionless rheological | Definition | Relevance | Values
---|---|---|---
parameter |  |  |
$Cu$ | $Cu=\lambda_{t}\dot{\gamma_{c}}$ | Ratio of characteristic strain rate | $0\leq Cu<\infty$
|  | and critical strain rate |
$n$ | Power law index | Slope of shear thinning region | $0<n<1$
$\beta$ | $\beta=\frac{\eta_{\infty}}{\eta_{0}}$ | Ratio between infinite shear viscosity | $\beta\ll 1$
|  | and zero shear viscosity |
$\epsilon$ | $1-\beta$ | Relative viscosity drop | $\epsilon\sim 1$

####  2.2.3 Perturbation expansion for constitutive equation

The rheological model under consideration (Eq. (4)) is nonlinear in
$\dot{\gamma},Cu$ and $n$. In anticipation of the subsequent
microhydrodynamical study, we will express the viscosity in the following form
[5]:

$\eta=\mu+\delta A{[\dot{\gamma}]},$ (7)

where $\mu$ is a (constant) Newtonian viscosity, $A$ is a function of the
shear rate, and $\delta$ is a small parameter such that when $\delta\to 0$,
the Newtonian limit is recovered. We will examine two different perturbative
limits.

##### Small Carreau number:

When $Cu^{2}\dot{\gamma}^{2}\ll 1$, one can Taylor expand the viscosity around
the zero shear rate plateau. The small parameter is $\delta=Cu^{2}$ in the
perturbation expansion. We find $\mu=1$ and
$A[\dot{\gamma}]=\frac{1}{2}\epsilon(n-1)\dot{\gamma}^{2}$, which yields:

$\displaystyle\eta=1+\frac{1}{2}\epsilon(n-1)Cu^{2}\dot{\gamma}^{2}+O(Cu^{4}),$
(8)

##### Large Carreau number:

When $Cu^{2}\dot{\gamma}^{2}\gg 1$, one can Taylor expand the viscosity around
its infinite shear rate plateau. Here, the small parameter is
$\delta=Cu^{(n-1)}$ in the perturbation expansion. We find $\mu=\beta$ and
$A[\dot{\gamma}]=\epsilon\dot{\gamma}^{n-1}$, which yields:

$\eta=\beta+Cu^{(n-1)}\epsilon\dot{\gamma}^{n-1}+O(Cu^{2(n-1)})$ (9)

In Fig. 3, we compare the full rheological equation (Eq. (4)) to the
asymptotic limits discussed above. Overall, we see that the small Carreau
number limit does a reasonable job capturing the rheology for $Cu^{2}<0.1$,
while the large Carreau number limit does a reasonable job for $Cu>10$.

Figure 3: Viscosity of Carreau fluid for $n=0.3$ and $\epsilon=0.9$. The full
rheological equation (Eq. (4)) is compared against the perturbation results in
the small Carreau number limit (Eq. (8)) and the large Carreau number limit
(Eq. (9)).

##  3 Theory

###  3.1 Overview of steps and non-dimensionalization

The aim of this work is to analyse the microhydrodynamics of spheroids in
shear thinning fluids. We will segregate our analysis into two parts. First,
we will calculate the background flow field of the shear thinning fluid
without the particle. This analysis will be followed by an investigation of a
spheroid in this flow field. The reciprocal theorem will be used to obtain
corrections to the rigid body motion, in the limits of small and large Carreau
numbers.

From here on out, we will write results in non-dimensional form. All distances
will be scaled by the channel half height $h$. The viscosity will be scaled by
the zero-shear rate viscosity $\eta_{0}$, while the stresses will be scaled by
$\tau_{c}=\frac{h\Delta P}{2l}$. The shear rate will be scaled by
$\dot{\gamma}_{c}=\tau_{c}/\eta_{0}=\frac{h\Delta P}{2\eta_{0}l}$, the time
will be scaled by $\dot{\gamma}_{c}^{-1}$, while the velocities will be scaled
by $V_{c}=\dot{\gamma}_{c}h$. The Carreau number for this geometry will be
defined as:

$Cu=\frac{\lambda_{t}h\Delta P}{2\eta_{0}l}$ (10)

###  3.2 Background flow field

In this subsection, we will use the coordinate system in Figure 1 (a) and
solve for the flow field in the channel. The quantities in this subsection
will have a superscript “$\infty$” indicating that they are solved in the
absence of a particle.

When the flow is steady, inertialess, and incompressible, the momentum and
continuity equations state that $\frac{\partial\tau_{ij}^{\infty}}{\partial
x_{j}}=0$ and $\frac{\partial u_{i}^{\infty}}{\partial x_{i}}=0$ (Einstein
convention assumed). The flow is unidirectional (i.e.,
$u_{x}^{\infty}=u_{x}^{\infty}(y)$ only), which simplify the equations
considerably. The differential equation for the flow field is:

$\frac{\partial\tau_{xy}^{\infty}}{\partial y}=-2$ (11)

where $\tau_{xy}^{\infty}$ is given by the Carreau model (Eqns. (2) and (4))
with local shear rate
$\dot{\gamma}^{\infty}=\frac{\partial{u_{x}^{\infty}}}{\partial y}$. The
differential equation is subject to the boundary conditions $u^{\infty}_{x}=0$
at the top wall ($y=h$) and $\frac{\partial u_{x}^{\infty}}{\partial y}=0$ at
the center plane ($y=0$).

Below state the results for the velocity field in terms of the non-dimensional
quantities of the Carreau model (Table 1). These results were obtained by
performing a regular perturbation expansion in $Cu^{2}$ for small Carreau
number, and $Cu^{n-1}$ for large Carreau number. For small Carreau number
($Cu\ll 1$), the velocity field is:

$\displaystyle
u_{x}^{\infty}=u_{x}^{\infty,(0)}+Cu^{2}u_{x}^{\infty,(1)}+O(Cu^{4})$ (12a)
$\displaystyle u_{x}^{\infty,(0)}=1–y^{2};\qquad
u_{x}^{\infty,(1)}=\epsilon(1-n)(1–y^{4})$ (12b)

while for large Carreau number ($Cu\gg 1$) the velocity field is:

$\displaystyle
u_{x}^{\infty}=u_{x}^{\infty,(0)}+Cu^{n-1}u_{x}^{\infty,(1)}+O(Cu^{2(n-1)})$
(13a) $\displaystyle u_{x}^{\infty,(0)}=\frac{1–y^{2}}{\beta};\qquad
u_{x}^{\infty,(1)}=-\frac{2^{n}\epsilon}{\beta^{n+1}(n+1)}(1–|y|^{n+1})$ (13b)

For both cases, the shear stress is $\tau_{xy}=-2y$. Fig. 4 plots the flow
field in the two perturbative limits. The small Carreau approximation captures
the deviation from the first Newtonian plateau, where the normalized viscosity
$\mu=1$. A perturbation around that plateau gives rise to a reduced effective
viscosity, ultimately resulting in a higher flow rate. Compare this with the
large Carreau number approximation, which captures deviations around the
second Newtonian plateau with a normalized viscosity $\mu=\beta.$ A
perturbation around this plateau gives rise to a higher viscosity, which
reduces the flow rate.

(a)

(b)

Figure 4: Pressure driven flow of shear thinning fluids in limits of (a)
$Cu\ll 1$ (Eq. (12)) and (b) $Cu\gg 1$ (Eq. (13)). The rheological parameters
are $\epsilon=0.9$ and $n=0.3$.

###  3.3 Microhydrodynamics

When the spheroid’s size is much smaller than the channel height (i.e., $R\ll
h$), one can determine the motion of the spheroid by assuming it is in an
unbound fluid with a background velocity given by Eqs. (12) and (13). Here, we
will use the coordinate system in Figure 1(b) to compute the particle’s rigid
body motion. The coordinate $x_{i}$ represents the position from the
particle’s center of mass, while $u_{i}^{\infty}$ is the background flow from
the previous section, transformed to the particle coordinate system. We will
use the reciprocal theorem to approximate expressions for the rigid body
motion in the limit of small and large Carreau numbers.

####  3.3.1 Problem setup, reciprocal theorem, and perturbation expansion

Suppose we have an inertialess, incompressible fluid that satisfies the
Carreau constitutive relationship Eqs. (2). This stress tensor can be
decomposed into two components:

$\displaystyle\tau_{ij}$ $\displaystyle=\tau_{ij}^{N}+\tau_{ij}^{ex}$ (14a)
$\displaystyle\tau_{ij}^{N}$ $\displaystyle=\mu\dot{\gamma}_{ij}–p\delta_{ij}$
(14b) $\displaystyle\tau_{ij}^{ex}$
$\displaystyle=(\eta-\mu)\dot{\gamma}_{ij}$ (14c)

where $\tau_{ij}^{N}$ is a Newtonian stress tensor with constant viscosity
$\mu$, while $\tau_{ij}^{ex}$ is an extra stress tensor that depends on the
shear-thinning viscosity $\eta$ described in Eq. (4).

When we place a spheroid in such a fluid with far field velocity
$u_{i}^{\infty}$, the continuity and momentum equations outside the particle
can be written as:

$\frac{\partial u_{i}}{\partial
x_{i}}=0;\qquad\frac{\partial\tau_{ij}^{N}}{\partial x_{j}}+b_{i}=0$ (15)

where

$b_{i}=\frac{\partial}{\partial x_{j}}\left(\tau_{ij}^{ex}\right)$ (16)

is an effective body force on a Newtonian fluid. Far away from the particle
($|x_{i}|\rightarrow\infty$), the velocity and stress fields are
$u_{i}^{\infty}$, $\tau_{ij}^{N,\infty}$, and $\tau_{ij}^{ex,\infty}$, which
also satisfy the above relationship (15) throughout the entire domain. On the
surface $S_{p}$ of the particle, the velocity is rigid body motion:

$u_{i}=U_{i}+\epsilon_{ijk}\Omega_{j}x_{k}\qquad\text{on }\qquad x_{i}\in
S_{p}$ (17)

where $U_{i}$ and $\Omega_{i}$ are the translational and rotational velocities
of the particle. The force and torque balance also need to be satisfied. The
external force and torque on the particle are:

$\displaystyle F_{i}^{ext}$ $\displaystyle=\int_{S_{p}}\tau_{ij}n_{j}dS$ (18a)
$\displaystyle T_{i}^{ext}$
$\displaystyle=\int_{S_{p}}\epsilon_{ijk}x_{j}\tau_{km}n_{m}dS$ (18b)

where the normal vector $n_{i}$ points inward to the particle. In this
problem, $F_{i}^{ext}=T_{i}^{ext}=0$. The force and torque on the particle are
also zero from the far-field stress
$\tau_{ij}^{\infty}=\tau_{ij}^{N,\infty}+\tau_{ij}^{ex,\infty}$.

We note that the above problem is formulated as the Stokes flow around a
particle with a far field velocity $u_{i}^{\infty}$ and a fluid body force
$b_{i}$. In this situation, one can employ the reciprical theorem to obtain an
expression for the rigid body motion $(U_{i},\Omega_{i})$. This theorem has a
celebrated history in the Stokes flow community – see classical texts [4, 2].
Below, we provide the final formula for $(U_{i},\Omega_{i})$ – a detailed
derivation can be found in Appendix A.

The expression for the particle’s rigid body motion satisfies the following
relationship:

$\begin{bmatrix}R^{FU}_{ij}&R^{F\Omega}_{ij}\\\
R^{TU}_{ij}&R^{T\Omega}_{ij}\end{bmatrix}\begin{bmatrix}U_{j}\\\
\Omega_{j}\end{bmatrix}=\begin{bmatrix}F_{i}^{eff}\\\
T_{i}^{eff}\end{bmatrix}$ (19)

In the above equation, the quantities $R_{ij}^{FU}$, $R_{ij}^{F\Omega}$,
$R_{ij}^{TU}$ and $R_{ij}^{T\Omega}$ are the Stokes-flow resistance tensors
for the particle with Newtonian viscosity $\mu$. This formula states that the
the translation and rotation of the particle is equal to its motion in Stokes
flow with an effective force and torque:

$\displaystyle F_{i}^{eff}=F_{i}^{ext}+F_{i}^{flow}+F_{i}^{poly}$ (20a)
$\displaystyle T_{i}^{eff}=T_{i}^{ext}+T_{i}^{flow}+T_{i}^{poly}$ (20b)

The first terms in the above expression correspond to the external force and
torque ($F_{i}^{ext}$, $T_{i}^{ext}$) on the particle, which are zero for this
problem. The second terms correspond to the effective force and torque from
the external flow $u_{i}^{\infty}$. These quantities are:

$\displaystyle F_{k}^{flow}$
$\displaystyle=\int_{S_{p}}u_{i}^{\infty}\Sigma_{ijk}^{trans}n_{j}dS$ (21a)
$\displaystyle T_{k}^{flow}$
$\displaystyle=\int_{S_{p}}u_{i}^{\infty}\Sigma_{ijk}^{rot}n_{j}dS$ (21b)

In the above formula, the quantities $\Sigma_{ijk}^{trans}$ and
$\Sigma_{ijk}^{rot}$ are the the $ij$ components of the stress field on the
particle surface in Stokes flow due to unit translation or unit rotation in
the $k$ direction. The normal vector $n_{j}$ points into the particle, and the
expression is integrated over the particle surface $S_{p}$.

The last terms in Eq. (20) are the effective force and torque from the non-
Newtonian stress. These quantities are:

$\displaystyle F_{k}^{poly}$ $\displaystyle=-\int_{V}\frac{\partial
v_{ik}^{trans}}{\partial
x_{j}}\left(\tau_{ij}^{ex}-\tau_{ij}^{ex,\infty}\right)dV$ (22a)
$\displaystyle T_{k}^{poly}$ $\displaystyle=-\int_{V}\frac{\partial
v_{ik}^{rot}}{\partial
x_{j}}\left(\tau_{ij}^{ex}-\tau_{ij}^{ex,\infty}\right)dV$ (22b)

where the volume of integration is outside the particle. The quantities
$v_{ik}^{trans}$ and $v_{ik}^{rot}$ represent the Stokes velocity field
outside the particle in the $i$ direction due to unit translation or unit
rotation in the $k$ direction. The quantity $\tau_{ij}^{ex}$ is the extra
stress outside the particle, while $\tau_{ij}^{ex,\infty}$ is the extra stress
in the absence of the particle.

We will now perform a perturbation expansion to estimate the rigid body motion
$(U_{i},\Omega_{i})$ in the limits of small and large Carreau numbers. We will
expand the velocity fields as follows:

$\displaystyle u_{i}$ $\displaystyle=u_{i}^{(0)}+\delta
u_{i}^{(1)}+O(\delta^{2})$ (23a) $\displaystyle u_{i}^{\infty}$
$\displaystyle=u_{i}^{\infty,(0)}+\delta u_{i}^{\infty,(1)}+O(\delta^{2})$
(23b)

and similarly expand the viscosity of the fluid as:

$\eta=\mu+\delta A[\dot{\gamma}^{(0)}]+O(\delta^{2})$ (24)

In the above equations, $\delta$ is a small parameter, where $\delta=Cu^{2}$
for small Carreau number $(Cu\ll 1)$, and $\delta=Cu^{(n-1)}$ for large
Carreau number $(Cu\gg 1)$. For small Carreau number, the constant viscosity
is $\mu=1$, while the nonlinear function is
$A[\dot{\gamma}]=\frac{1}{2}\epsilon(n-1)\dot{\gamma}^{2}$. This situation
corresponds to a weak departure from the zero shear rate plateau. For large
Carreau number, the constant viscosity is $\mu=\beta$, while
$A[\dot{\gamma}]=\epsilon\dot{\gamma}^{n-1}$. This situation corresponds to a
weak departure from the infinite shear rate plateau.

We will now plug in the above expansions into Eqs. (19)-(22) to obtain the
particle’s translation and rotational velocity $(U_{i},\Omega_{i})$ to
$O(\delta)$. We specifically require the far-field velocity $u_{i}^{\infty}$
up to $O(\delta)$ for the flow force and torque (Eq. (21)). The expression for
$u_{i}^{\infty}$ is found in Eqs. (12) and (13), which needs to be transformed
to the particle coordinate system. We also require the extra stress tensor
$\tau_{ij}^{ex}-\tau_{ij}^{ex,\infty}$ up to $O(\delta)$ for the polymeric
force and torque (Eq. (22)). The expression for this quantity is:

$\tau_{ij}^{ex}-\tau_{ij}^{ex,\infty}=\delta
A[\dot{\gamma}^{(0)}]\dot{\gamma}_{ij}^{(0)}-\delta
A[\dot{\gamma}^{\infty,(0)}]\dot{\gamma}_{ij}^{\infty,(0)}+O(\delta^{2})$ (25)

where as mentioned before,
$A[\dot{\gamma}]=\frac{1}{2}\epsilon(n-1)\dot{\gamma}^{2}$ in the small
Carreau number limit ($Cu\ll 1)$, while
$A[\dot{\gamma}]=\epsilon\dot{\gamma}^{n-1}$ in the large Carreau number limit
$(Cu\gg 1)$. The above expression depends on the $O(1)$ shear rate tensors
$\dot{\gamma}^{(0),\infty}_{ij}$ and $\dot{\gamma}^{(0)}_{ij}$, which are
spatial derivatives of the $O(1)$ velocity fields. The $O(1)$ velocity field
$u_{i}^{(0)}$ is the Stokes flow around the particle with a far-field velocity
$u_{i}^{\infty,(0)}$. For the specific case where the particle is an ellipsoid
and the far-field velocity is parabolic (e.g., Eqs. (12) and (13) at $O(1)$),
the solution to $u_{i}^{(0)}$ is known [47, 48].

####  3.3.2 Formulae for ellipsoids

Eqs. (19)-(22) describe the rigid body motion of a particle in a non-Newtonian
fluid. For a specified velocity field $u_{i}^{\infty}$ and extra stress tensor
$\tau_{ij}^{ex}$, one can compute the particle’s translational and rotational
velocity $(U_{i},\Omega_{i})$. Performing this calculation requires knowledge
of the particle’s Stokes flow behavior. Specifically, one needs to know the
velocity fields $(v_{ik}^{trans},v_{ik}^{rot})$ and stress fields
$(\Sigma_{ijk}^{trans},\Sigma_{ijk}^{rot})$ outside the particle arising from
rigid body motion, and one needs to know the resistance tensors
($R_{ij}^{FU}$, $R_{ij}^{F\Omega}$, $R_{ij}^{TU},R_{ij}^{T\Omega}$). These
quantities are well-known for an ellipsoid.

Let us consider an ellipsoid with semi-axes $(a,b,c)=(a_{1},a_{2},a_{3})$, and
let us choose a particle coordinate system that aligns with these axes. In
this situation, Eqs. (19)-(22) simplify considerably for the particle’s rigid
body motion. For the case when the external force and torque are zero
($F_{i}^{ext}=T_{i}^{ext}=0$), we obtain:

$\displaystyle U_{i}$
$\displaystyle=U_{i}^{Faxen}+\frac{1}{R_{ii}^{FU}}F_{i}^{poly}\qquad\text{(no
summation over i)}$ (26a) $\displaystyle\Omega_{i}$
$\displaystyle=\Omega_{i}^{Faxen}+\frac{1}{R_{ii}^{T\Omega}}T_{i}^{poly}\qquad\text{(no
summation over i)}$ (26b)

where $U_{i}^{Faxen}$ and $\Omega_{i}^{Faxen}$ are the translational and
rotational velocities from Faxen’s formula for an ellipsoid. These quantities
are:

$\displaystyle
U_{i}^{Faxen}=-\frac{1}{3V_{p}}\int_{S_{p}}u_{i}^{\infty}(n_{j}x_{j})dS$ (27a)
$\displaystyle\Omega_{i}^{Faxen}=-\frac{1}{V_{p}}P_{ij}\epsilon_{jrk}\int_{S_{p}}x_{r}u_{k}^{\infty}(n_{m}x_{m})dS$
(27b)

where $V_{p}=\frac{4\pi}{3}abc$ is the particle volume and $P_{ij}$ is a
diagonal tensor. The 11 component of this tensor is
$P_{11}=\frac{1}{a_{2}^{2}+a_{3}^{2}}$ with the other diagonal components
obtained from index cyling. If we Taylor expand the velocity field around the
center of mass $(x_{i}=0)$, one can get an approximate expression for the
Faxen’s velocities up to $O((R/h)^{2})$, which is within the current current
approximation of our model. This procedure yields:

$U_{i}^{Faxen}=u_{i}^{\infty}(0)+\frac{1}{6}\sum_{k=1}^{3}a_{k}^{2}\frac{\partial
u_{i}^{\infty}}{\partial x_{k}\partial x_{k}}(0)+\dots$ (28)

where the far field velocity $u_{i}^{\infty}$ and its derivatives are
evaluated at the center of mass. For the rotational velocity, one obtains:

$\Omega_{1}^{Faxen}=\frac{1}{2}\omega_{1}^{\infty}(0)+\frac{a_{2}^{2}-a_{3}^{2}}{a_{2}^{2}+a_{3}^{2}}E_{23}^{\infty}(0)+\dots$
(29)

where $\omega_{i}^{\infty}=\epsilon_{ijk}\frac{\partial
u_{k}^{\infty}}{\partial x_{j}}$ is the vorticity evaluated at the center of
mass, and $E_{ij}^{\infty}=\frac{1}{2}\dot{\gamma}_{ij}^{\infty}$ is the rate
of strain tensor evaluated at the center of mass. The other components for the
rotational velocity are obtained by index cycling.

In Eq. (26), the resistance tensors $R_{ii}^{FU}$ and $R_{ii}^{T\Omega}$ take
the following form:

$\displaystyle R_{11}^{FU}=16\pi\mu abc\frac{1}{\chi_{0}+\alpha_{1}a_{1}^{2}}$
(30a) $\displaystyle R_{11}^{T\Omega}=\frac{16\pi\mu
abc}{3}\frac{a_{2}^{2}+a_{3}^{2}}{a_{2}^{2}\alpha_{2}+a_{3}^{2}\alpha_{3}}$
(30b)

In the above equations, $\chi_{0}$ and $(\alpha_{1},\alpha_{2},\alpha_{3})$
are elliptic integrals given in publication [48]. If one wants to obtain the
other components of the diagonal tensor, one changes the indices for $a_{i}$
and $\alpha_{i}$ appropriately. For example, for for $R_{22}^{FU}$, one
changes the index to “2” while for $R_{22}^{T\Omega}$ one changes the indices
to “1” and “3”.

To evaluate the polymeric force $F_{i}^{poly}$ and torque $T_{i}^{poly}$ on
the particle, one uses Eq. (22), which uses the Stokes velocity fields
$v_{ik}^{trans}$ and $v_{ik}^{rot}$ outside the particle due to unit
translation and unit rotation. Full expressions for $v_{ik}^{trans}$ and
$v_{ik}^{rot}$ are given in Appendix B. This formula also requires knowledge
of the extra stress tensor $\tau_{ij}^{ex}-\tau_{ij}^{ex,\infty}$ (Eq. (25)),
which involves the $O(1)$ velocity field around the particle. The solution to
the $O(1)$ velocity field around an ellipsoid is calculated in papers [47, 48]
(see also [16]). We will refer readers to these works for the full
expressions.

##  4 Code development and verification

Given a particle’s initial position $(x_{0},y_{0},z_{0})$, and orientation
($\theta_{0},\phi_{0}$) in the channel (Fig. 1 a), we track its position and
orientation over a time as follows:

  1. 1.

Transform from lab coordinates to particle coordinates: We transform from lab
coordinates (Fig.1 (a)) to particle coordinates (Fig. 1 (b)), where the origin
is now at the spheroid’s center of mass and the axes are aligned with the
spheroid’s semi-axes. We transform the channel velocity field $u_{i}^{\infty}$
(Eqs. (12) and (13)) to this coordinate system as well. This velocity field is
evaluated up to $O(\delta)$, where $\delta=Cu^{2}$ for the small Carreau
number limit and $\delta=Cu^{n-1}$ for the large Carreau number limit.

  2. 2.

Evaluate Faxen velocities: Using the expression for $u_{i}^{\infty}$, we
evaluate the Faxen velocities $U_{i}^{Faxen}$ and $\Omega_{i}^{Faxen}$ using
Eqns. (28) and (29).

  3. 3.

Calculate polymer force and torque: We calculate the polymer force and torque
on the spheroid using Eq. (22). For each component of the polymer force and
torque, we numerically evaluate the volume integral in ellipsoidal coordinates
using Gaussian quadrature. We typically use $15$ quadrature points in the
radial direction, $10$ quadrature points in the latitude direction, and $35$
quadrature points in the longitude direction. In the integrand, we use the
expression for the extra stress $\tau_{ij}^{ex}-\tau_{ij}^{ex,\infty}$
described in Eq. (25). Note that the integrand also involves the quantities
from the auxiliary problem $(v_{ik}^{trans},v_{ik}^{rot})$, which has been
solved and tabulated separately beforehand.

  4. 4.

Calculate rigid body motion: Using the Faxen velocities and the polymer
force/torque, we calculate the rigid body motion of the particle in the
particle coordinate frame (Eqs. (26) and (30)).

  5. 5.

Transform back into lab coordinates and update particle position and
orientation: We transform the rigid body velocities $(U_{i},\Omega_{i})$ back
to the lab coordinates. We then update the particle’s position and
orientation. The equation governing the particle’s center of mass is:

$\frac{dx_{i}^{cm}}{dt}=U_{i}$ (31)

The equation governing the time evolution of the projection vector is
$\frac{dp_{i}}{dt}=\epsilon_{ijk}\Omega_{j}p_{k}$. If written in terms of the
orientation angles $(\theta,\phi)$ of the spheroid, these equations become:

$\begin{gathered}\frac{d\theta}{dt}=\Omega_{y}\cos\phi-\Omega_{x}\sin\phi\\\
\frac{d\phi}{dt}=\Omega_{z}-\Omega_{x}\cos\phi\cot\theta-\Omega_{y}\sin\phi\cot\theta,\end{gathered}$
(32)

We update the particle’s position and orientation angles using the above
equations (31) and (32). We perform forward Euler time-stepping, with the time
step chosen between $10^{-3}<\Delta t<10^{-2}$ in dimensionless units.

  6. 6.

Repeat the above steps: Repeat steps 1-5 until we have simulated the particle
motion over a sufficient period of time.

We developed a MATLAB code to execute this algorithm. The code was initially
developed and used for the analysis of an ellipsoid in a second order
viscoelastic fluid in a previous publication [16]. For this paper, we modified
the code to account for shear thinning effects in both the low Carreau and
high Carreau number limits. To test the code, we simulated the sedimentation
of a falling sphere in a shear thinning fluid in the low Carreau number limit.
Fig. 5 show the results of our simulation and theoretical results provided by
[24]. The exact match between the results serves to verify our code.

Figure 5: Sedimentation velocity of a sphere in a shear thinning fluid. The
$y$-axis shows the fractional increase of the sedimentation speed $U$ compared
to the Newtonian value $U_{N}$ at zero Carreau number ($Cu=0$). The solid
lines are results from the theory of [24], while the dots are the results of
the reciprocal theorem based numerical simulation carried out in this paper.

##  5 Results and discussion

In our discussion, we will stick to the orientational dynamics of spheroids
and forego a discussion of translational motion. We find that for
translational motion, no novel phenomena are revealed like cross stream
migration or particle lift [24, 25]. We will first discuss the results in the
small Carreau number limit, followed by the results in the large Carreau
number limit.

Unless otherwise noted, we will examine a particle with equivalent particle
radius $R=0.1$ when rendered dimensionless by the channel height $h$. The
initial position of the particle is at $[x_{0},y_{0},z_{0}]=[0,-0.5,0]$.

###  5.1 Orientation dynamics – small Carreau number limit ($Cu\ll 1$)

####  5.1.1 General observations – effect of $A_{R}$, $Cu$, and $n$:

(a)

(b)

(c)

(d)

Figure 6: Effect of shape parameter ($A_{R}$) on the orientation trajectories
of prolate and oblate spheroids in the small Carreau number limit $(Cu\ll 1)$.
In all these plots, the initial orientation is
$\theta_{0}=\pi/4,\phi_{0}=\pi/4$. The plots on the left show prolate
spheroids ($A_{R}>1$) and those on the right show oblate spheroids
($A_{R}<1$). The solid curves denote shear thinning fluids with $Cu^{2}=0.1$
and $n=0.3$, while the dashed curves denote Newtonian fluids ($Cu^{2}=0$).

Fig. 6 shows the influence of the shape parameter ($A_{R}$) on the
orientational dynamics of spheroids in Newtonian and shear thinning fluids.
First and foremost, we remark that the further a particle deviates from a
spherical shape, the larger is the time period for its tumbling. For a prolate
particle ($A_{R}>1$), the time period increases as $A_{R}$ increases. For an
oblate particle $(A_{R}<1)$, the time period increases as $A_{R}$ decreases.
These trends are clearly an artefact of Newtonian flow, where the time period
scales as $T\sim(A_{R}+1/A_{R})$. Next, we observe that the shear thinning
reduces the tumbling time period compared to a Newtonian fluid —- see the
solid (shear thinning) and dashed (Newtonian) curves in Fig. 6. Additionally,
shear thinning reduces the amplitude of oscillations in $\theta$ for the
particles. These trends are contrary to the trends observed for prolate
spheroids in linear flows of shear thinning fluids, as shown in the $Cu^{2}\ll
1$ limit in the work of [25]. In their case, shear thinning increases the time
period of revolution for prolate spheroids, and increases the oscillation
amplitude for $\theta$. This dichotomy in the trends is explained below.

For pressure driven flows analysed here, the deviatoric stress field balances
the pressure gradient across the channel imposed on its ends. Since the
imposed pressure drop is independent of the shear thinning correction, the
stress induced in the flow does not change due to shear thinning. Since the
stress field induced in the fluid is the same as that in the Newtonian case, a
shear thinning fluid, being less viscous, has to deform more than a Newtonian
fluid to sustain the same stress and pressure drop. Because the time period
$T$ and the oscillation amplitude for $\theta$ scale as the inverse of strain
rate/deformation, both these quantities decrease with an increase in shear
thinning as shown in Fig. 6.

On the other hand, the linear (shear) flows studied by [25] are kinematically
controlled, wherein the velocity field of the background flow remains
unchanged between Newtonian and shear thinning rheologies. Therefore, the
stress field induced in a shear thinning flow field is weaker than that of a
Newtonian field. A weaker stress field in turn imposes a diminished
hydrodynamic torque on the particle, slowing down its rotation and increasing
its oscillation amplitude $\theta$. The discussion of tumbling behavior of
spheroids in shear flows and their comparison with the pressure driven flows
is further elaborated upon in Sec. 6.

Fig. 7 plots the tumbling time period for prolate and oblate spheroids for
different values of Carreau number ($Cu$), power law index ($n$) and shape
parameter ($A_{R}$). A couple of trends can be noticed here. First, we clearly
see that increased shear thinning (i.e., an increase in $Cu$ or a decrease in
$n$) reduces the time period of tumbling, and this effect is amplified the
further the particle deviates from a sphere. Even in the case of linear shear
flows in [25], the effects of shear thinning are similarly enhanced at higher
aspect ratios for prolate particles. The straightforward reason is that as the
aspect ratio ($A_{R}$) of a prolate particle increases (for the same volume),
the projection of the long axis along the velocity gradient $y$-axis (given by
$a\sin{\theta}\sin{\phi}$) also becomes larger, and the particle is exposed
more strongly to the shear thinning tendency of the flow. A similar argument
holds for oblate particles as well. Interestingly, we see that time period is
exactly equal when comparing the $A_{R}=3,4,5$ prolate particles to the
$A_{R}=\frac{1}{3},\frac{1}{4},\frac{1}{5}$ oblate particles. This observation
suggests that the just like Newtonian case, the time period ${T}$ is a
function of $(A_{R}+1/A_{R})$.

The last point we would like to illustrate is that the decrease in the
tumbling time period cannot be wholly explained by the reduction in the
effective shear rate around the particle. In a Newtonian flow, the time period
of tumbling is given by the classical Jeffrey formula
$T=2\pi(A_{R}+1/A_{R})\dot{\gamma}_{\text{loc}}^{-1}$, where the
$\dot{\gamma}_{\text{loc}}$ is the local strain rate at the particle’s center
of mass $y=y_{0}$. A naïve approach to estimating the tumbling time period in
a shear thinning fluid would be to use the same formula, noting that
$\dot{\gamma}_{loc}$ increases in the channel as $Cu$ increases or $n$
decreases. Fig. 7 shows the time period estimated using this approach, plotted
as dotted lines. Overall, while we see this approach accounts for some of the
decrease in the particle’s time period, it does not match the simulations
well. This indicates that one cannot use simple modifications of Jeffrey’s
theories to predict the time period of spheroids in shear thinning flows.

(a)

(b)

Figure 7: Tumbling time period for both prolate and oblate particles in the
$Cu\ll 1$ limit for different values of (a) shape parameter $A_{R}$ and (b)
power law index $n$. The initial orientation is
$\theta_{0}=\pi/4,\phi_{0}=\pi/4$ for both the plots. The solid lines show the
results from the reciprocal theorem based numerical estimation while the
dotted lines show the result using Jeffrey’s formula using the local shear
rate – i.e., $T=2\pi(A_{R}+1/A_{R})\dot{\gamma}_{\text{loc}}^{-1}$, where
$\dot{\gamma}_{\text{loc}}=\frac{\partial}{\partial
y}\left[u_{x}^{\infty,(0)}+Cu^{2}u_{x}^{\infty,(1)}\right]$

####  5.1.2 Effect of initial orientation ($\theta_{0},\phi_{0}$):

(a)

(b)

(c)

(d)

Figure 8: Effect of initial orientation ($\theta_{0},\phi_{0}$) on the
orientation trajectories of prolate and oblate spheroids in the small Carreau
number limit $(Cu\ll 1)$. The plots on the left show prolate spheroids with
$A_{R}=5$ and those on the right show oblate spheroids with $A_{R}=0.2$. The
solid lines denote shear thinning flows with power law index $n=0.3$ and
$Cu^{2}=0.1$, while the dashed lines denote Newtonian fluids ($Cu^{2}=0$).

Next, we analyse the effect of initial orientation ($\theta_{0},\phi_{0}$) on
the spheroid’s trajectory in both Newtonian and shear thinning flows. There
are six pairs of initial orientations we study: $(\theta_{0},\phi_{0})=$
($(\pi/4,0)$, $(\pi/4,\pi/4)$, $(\pi/4,\pi/2)$, $(\pi/3,0)$, $(\pi/3,\pi/4)$,
$(\pi/3,\pi/2)$). Both prolate and oblate particles are analysed.

First we observe from Fig. 8 (a,b) that the evolution of the azimuth angle
$\phi$ over time is independent of the initial orientation $\theta_{0}$ –in
other words, curves for both $\theta_{0}=\pi/4$ and $\theta_{0}=\pi/3$
coincide for the same initial angle $\phi_{0}$. For the prolate particle, the
spheroid spends more time in $\phi=m\pi$ configuration ($m=0,1$), while the
oblate particle spends more time in the $\phi=(m+\frac{1}{2})\pi$
configuration ($m=0,1$). The regions with $\phi=(m+\frac{1}{2})\pi$ (prolate
particles) and $\phi=m\pi$ (oblate particles) show the largest deviation
between Newtonian and shear thinning trajectories. We can understand these
trends by noting that the rotation of the spheroid is actuated by the
hydrodynamic torque from the shear stresses. This torque arises from a
traction in the $x$-direction and a moment arm (projection of particle length)
along the $y$-axis. Thus, configurations with the smallest moment arm (e.g.,
$\phi=m\pi$ for prolate particle, $\phi=(m+\frac{1}{2})\pi$ for oblate
particle) have the smallest rotation rate and longest residence time, while
configurations with the longest moment arm (e.g., $\phi=(m+\frac{1}{2})\pi$
for prolate particle, and $\phi=m\pi$ for the oblate particle) have the
largest rotation rate and hence largest influence of shear thinning.

At this juncture, we would like to compare the evolution of $\phi$ in this
flow with spheroids in other kinds of flow. For prolate spheroids in pure
shear flow of shear thinning fluids, the trend is same as observed here
wherein the prolate spheroid rotates the fastest at $\phi=(m+\frac{1}{2})\pi$
and slowest in the $\phi=m\pi$ configuration (see Fig.4 (b) and Fig.4 (c) in
the [25]). On the other hand, for the case of spheroids in quadratic flow of
nonlinear viscoelastic fluids, the reverse is true. For prolate (oblate)
spheroids in such flows, the maxima (minima) of $\dot{\phi}$ occurs at
$\phi=0$ or $\phi=\pi$, while the minima (maxima) occurs at $\phi=\pi/2$ (see
Fig. $6$ and the relevant discussion in Sec. IV c in [16]). The reason behind
this dichotomy is explained by considering the fact that in nonlinear
viscoelastic flows, the hydrodynamic torque acting on the spheroid is induced
by normal stresses, as opposed to shear stresses in shear thinning flows.
These normal stresses give rise to a traction in the $y$-direction and a
moment arm in the $x$-direction. Consequently, when the spheroid is oriented
such that its longest axis is aligned with the $x$-axis, the hydrodynamic
torque acting on the spheroid is the largest, resulting in higher rotation
rate in this configuration.

We next analyse the evolution of $\theta$ in time for both Newtonian and shear
thinning flows (Fig. 8 (c,d)). From the classical treatment of a Jefferey
orbit in a Newtonian fluid, we understand that $\dot{\theta}$ varies as
$\frac{A_{R}^{2}-1}{A_{R}^{2}+1}\sin{2\theta}\sin{2\phi}$. We also observe
similar features of this equation in shear thinning fluids. For example, we
see $\dot{\theta}$ change between prolate and oblate particles of the same
initial orientation due to the sign change in
$\frac{A_{R}^{2}-1}{A_{R}^{2}+1}$ (compare the same color curves in the Fig. 8
(c,d)). We also see that a spheroid initially released at $\theta_{0}=\pi/2$
or $\theta_{0}=0$ will have $\dot{\theta}=0$, and therefore remain at the
initial angle. Overall, we see shear thinning reduces the period and amplitude
of $\theta$; however, the influence of the initial orientation on these
quantities is seen to be less dominant than the other effects discussed
earlier ($A_{R}$, $Cu$, and $n$) . Therefore, we do not comment on this effect
further in the manuscript.

###  5.2 Orientation dynamics – large Carreau number limit ($Cu\gg 1$)

We now analyse the rotational dynamics of a spheroid in the large Carreau
number limit ($Cu\gg 1)$, which captures deviations in the fluid viscosity
from the lower Newtonian plateau $\eta=\beta$. There are a couple of
differences between the mathematical behavior in the $Cu\ll 1$ limit and the
$Cu\gg 1$ limit. First, since the $Cu\gg 1$ limit captures the deviation from
the second Newtonian plateau, any perturbation from this plateau will increase
the effective viscosity of the fluid. Secondly, for the large Carreau number
limit, the second Newtonian plateau is recovered solely in the limit
$Cu\to\infty$, but not in the limit $n\to 1$. This is in contrast to the small
Carreau number limit where the upper pleateau is recovered by $n\to 1$ or
$Cu\to 0$ or both.

####  5.2.1 General observations – effect of $A_{R}$, $Cu$, and $n$

(a)

(b)

(c)

(d)

Figure 9: Effect of shape parameter ($A_{R}$) on the orientation trajectories
of prolate and oblate spheroids in the large Carreau number limit $(Cu\gg 1)$.
The plots on the left show prolate spheroids and those on the right show
oblate spheroids. The solid lines denote shear thinning flows with $Cu=1000$
and $n=0.3$, while the dashed lines denote the Newtonian plateau
($Cu\rightarrow\infty$) with dimensionless viscosity $\beta=0.1$. For all the
plots, the initial orientation is $\theta_{0}=\pi/4,\phi_{0}=\pi/4$.

Fig. 9 plots the angle trajectories $(\phi(t),\theta(t))$ for oblate and
prolate particles at the second Newtonian pleateau ($Cu\rightarrow\infty$
limit) and the shear-thinning region ($Cu\gg 1$ but finite). Different values
of shape parameter ($A_{R}$) are considered. A few key observations can be
made. First, the figure shows that the time period and the amplitude for
$\theta$ are smaller in the Newtonian plateau regime ($Cu\rightarrow\infty$)
than the shear thinning regime ($Cu$ finite). The reason behind these trends
is that the shear thinning regime has a larger viscosity than the Newtonian
plateau, which leads to lower strain rates in the fluid, and hence larger time
periods for tumbling and larger oscillation amplitude for $\theta$ as
discussed previously. The figure also shows that the further the particle
shape deviates from a sphere ($A_{R}>1$ for prolates and $A_{R}<1$ for
oblates), the difference in dynamics amplifies between the Newtonian plateau
case and the shear thinning case. Lastly, orientation-wise, the angle $\phi$
behaves similarly in both the large Carreau and the small Carreau number
limits (compare Fig. 9 with Fig. 6). In both cases, the spheroid likes to
spend more time with the long axis aligned in the flow ($x$) direction, which
corresponds to $\phi=m\pi$ for prolate particles and $\phi=(m+\frac{1}{2})\pi$
for oblate particles. In these orientations, the difference in shear thinning
and Newtonian cases is minimum (see Fig. 9(a,b)), while the orientations with
the largest moment arm exhibit the largest differences between the shear
thinning and Newtonian cases.

Fig. 10 plots the tumbling time period of spheroids in the large Carreau
number regime for different values of the shape parameter ($A_{R}$), Carreau
number ($Cu$), and power-law index ($n$). These plots show that the time
period decreases to a constant value as the Carreau number increases, reaching
a Newtonian plateau independent of $Cu$ as $Cu\rightarrow\infty$. The changes
in the period however are much more modest than what was seen in the small
Carreau number limit (compare Fig 7 to Fig. 10). The reason for this
observation is that changes in the viscosity scale as $\epsilon(n-1)Cu^{2}$ in
the small Carreau number limit, while the changes scale as $\epsilon Cu^{n-1}$
in the large Carreau number limit, which exhibits a weaker variation with
respect to both $Cu$ and $n$. Similar to what was observed before, we find the
time period is the same for prolate and oblate particles of the same aspect
ratio, and the effects of shear thinning are enhanced the further the
particles deviate from a sphere. Lastly, in Fig. 10 we also plot the period
computed from the classic Jeffrey orbit equation
$T=2\pi(A_{R}+1/A_{R})\dot{\gamma}_{loc}^{-1}$, where $\dot{\gamma}_{loc}$ is
the local shear rate at the center of mass of the particle. One would expect
that $\dot{\gamma}_{loc}$ will change with Carreau number ($Cu$) and power-law
index ($n$), and thus explain some of the trends seen in the graph. We see
that while the equation does explain some of the trends, it does not overlap
with the simulation results, although the agreement appears to be better than
the small Carreau number limit in Fig 7. Again, these results suggest that one
cannot use simple ideas from Stokes flow to model the tumbling behavior of
particles in shear thinning fluids. The last point we would like to illustrate
is that the time period in Fig. 10 has different axes than plot in the small
Carreau number limit (Fig. 7). The tumbling period in the
$Cu\rightarrow\infty$ limit is $\beta$ times smaller than the $Cu=0$ limit,
due to the fact that the shear rate in the channel is $1/\beta$ times larger.
This will play a role in the trends discussed next section when we compare the
tumbling behavior of spheroids in pressure driven flows versus simple shear
flow.

(a)

(b)

Figure 10: Tumbling time period for both prolate and oblate spheroids in the
$Cu\gg 1$ limit for different values of (a) shape parameter $A_{R}$ and (b)
power law index $n$. The initial orientation is
($\theta_{0}=\pi/4,\phi_{0}=\pi/4$) and the rheological parameter $\beta=0.1$.
For (a) $n=0.3$, and for (b) $A_{R}=5$. The solid lines show the results from
the reciprocal theorem based numerical estimation while the dotted lines show
the result using Jeffrey’s formula using the local shear rate – i.e.,
$T=2\pi(A_{R}+1/A_{R})\dot{\gamma}_{\text{loc}}^{-1}$, where
$\dot{\gamma}_{\text{loc}}=\frac{\partial}{\partial
y}\left[u_{x}^{\infty,(0)}+Cu^{n-1}u_{x}^{\infty,(1)}\right]$.

##  6 Comparison with microhydrodynamics of spheroids in simple shear flow

In this section, we explore the motion of a spheroid in a Couette flow of a
shear thinning fluid. Such a flow is shown schematically in Fig. 11, where as
before, the coordinate axes are placed in the center of the channel. However,
here the flow is actuated by the top wall moving at a constant velocity
$2\mathcal{U}_{c}$. We non-dimensionalize all lengths by height $h$,
velocities by $\mathcal{U}_{c}$, shear rates by
$\dot{\gamma}_{c}=\mathcal{U}_{c}/h$, times by $\dot{\gamma}_{c}^{-1}$,
viscosities by $\eta_{0}$, and stresses by
$\tau_{c}=\eta_{0}\dot{\gamma}_{c}$.The velocity field in the absence of the
particle is given in dimensionless terms as:

Figure 11: Spheroid in a Couette flow of a shear thinning fluid. The apparatus
is same as that for pressure driven flow with the difference being that the
flow is actuated by the motion of the top wall.

${u}_{x}^{\infty}=1+y$ (33)

and the Carreau number is given by:

$Cu=\lambda_{t}\dot{\gamma}_{c}=\lambda_{t}\mathcal{U}_{c}/h$ (34)

We observe that the undisturbed velocity field is independent of viscosity,
and is therefore same for both Newtonian as well as for the perturbative shear
thinning cases ($\mathcal{O}(Cu^{2})$ and $\mathcal{O}(Cu^{n-1}$)).

The introduction of a spheroid into the linear background flow perpetuates a
disturbance field around the particle, in addition to the undisturbed
background flow. As before, to delineate the particle motion, we perturb the
disturbance field in two limits: small Carreau number (with $Cu^{2}$ as the
perturbation parameter) and large Carreau number (with $Cu^{n-1}$ as the
perturbation parameter). In both these perturbative solution schemes, the
leading order problem is Newtonian, whose solution has already been known in
literature and is available in standard textbooks ([49, Ch. 3]) For the shear
thinning correction, we appeal to the reciprocal theorem and solve for the
particle’s rigid body motion at $\mathcal{O}(Cu^{2})$ for small Carreau number
limit and at $\mathcal{O}(Cu^{n-1})$ for the large Carreau number limit.

The result of this analysis leads us to the tumbling period of a spheroid in
Couette flow of a shear thinning fluid in both the small Carreau number and
large Carreau number perturbative schemes. The results have been plotted in
Fig. 12. Here, we observe that for the small $Cu$ case in Fig 12(a) , the time
period increases with an increase in shear thinning. This result has also been
reported earlier by [25]. In the Newtonian limit of viscosity ($Cu\to 0$), the
time period also reduces to that given by the Jefferey orbit formula (see
dotted curve). More interesting for the current analysis is the result
pertaining to the large $Cu$ number perturbative limit, as plotted in
Fig.12(b). This is a novel result not published before in literature. Here, we
observe that as the fluid around the particle continues to shear thin, the
time period begins to decrease. At extremely large Carreau numbers, the time
period then decreases back to its Jefferey orbit formula. We note that in the
Newtonian limit of both these perturbative schemes ($Cu\to 0$ and
$Cu\to\infty$), the time period reduces to the same constant value.

(a)

(b)

Figure 12: Time period of revolution of prolate and oblate spheroids in linear
flows of shear thinning fluids in both (a) small $Cu$ and (b) large $Cu$
regimes. The solid curves show the shear thinning solution obtained using
reciprocal theorem based simulations, while the dashed curves show the
Newtonian solution using Jeffrey orbit formula. Rheological parameters include
$n=0.3$ and $\beta=0.1$ while the initial orientation is given by
$(\theta,\phi)=(\pi/4,\pi/4)$ Figure 13:  Cartoon of the time period of
revolution for spheroids in linear (top) and pressure driven (bottom) flows of
shear thinning fluids (figure not to scale). The solid red curves illustrate
known behavior from the perturbative analysis based on the reciprocal theorem,
while the dotted red curves are hypothesis. The dotted black lines show the
corresponding Newtonian plateaus. $T_{Newtonian2}=\beta\times T_{Newtonian1}$

In comparison to the pressure driven flow as discussed in the previous
sections, we notice the following differences in the microhydrodynamics of
spheroids in linear flows of shear thinning fluids:

  1. 1.

In pressure driven flows, the time period decreases monotonically with an
increase in Carreau number. On the other hand, in the case of linear flows,
the time period shows a non-monotonic behavior. In the low Carreau number
regime, the time period increases with increase in Carreau number; whilst on
the other hand, in the high Carreau number case, the time period decreases
with increases in Carreau number.

  2. 2.

In the pressure driven flows, the time period at the two Newtonian limits (low
and high Carreau) are different because the background flow and the Newtonian
viscosity at these two limits are different. In the linear flow, the time
period at the two Newtonian limits are the same, because the background flow
at these two limits are also the same.

These two differences are schematically shown in Fig. 13

##  7 Conclusion

In this paper, we analysed the motion of a spheroid in a pressure driven flow
of a shear thinning fluid. The shear thinning rheology is captured by the
Carreau model. We employ a perturbative approach in conjunction with the
reciprocal theorem to delineate the orientational dynamics of prolate and
oblate spheroids. First, an approximate version of the Carreau model is
rendered, in both the small and large Carreau number limits. Next, the
classical pressure driven (Poiseuille) flow equations are amended to take into
account the shear thinning corrections due to rheology. Finally, using the
reciprocal theorem and mobility relationships for ellipsoids in Newtonian
flow, the orientational kinematics of the problem is evaluated numerically.

We have the following conclusions to make:

  * •

Spheroids in shear thinning fluids tumble in a periodic fashion. In other
words, shear thinning does not resolve the degeneracy of Jeffrey’s orbits
found in Newtonian fluids. The degeneracy of Jefferey’s orbits is attributed
to the symmetry of the momentum equations for a Newtonian fluid. The shear
thinning rheology, even though nonlinear, still preserves the symmetry of the
momentum equations.

  * •

In the small Carreau number limit ($Cu\ll 1$), we find the spheroid’s tumbling
behavior to be very different in a pressure driven flow versus a linear shear
flow. In a pressure driven flow, shear thinning ($Cu>0$) gives rise to a
smaller tumbling time period compared to a Newtonian fluid ($Cu=0$). When the
particle is in a shear flow, shear thinning gives rise to the opposite trend –
i.e., a larger tumbling time period compared to its Newtonian counterpart
[25]. The reason for these opposing trends is that in the pressure driven flow
studied here, the flow is pressure controlled, while in a shear flow, the flow
is kinematically controlled. In the former case, shear thinning gives rise to
a larger shear rate on the particle, and hence a faster tumbling period. In
the latter case, shear thinning gives rise to a lower stress on the particle,
and hence a slower tumbling period.

  * •

We also observe a very different rotational behavior between spheroids in a
shear thinning fluid and a viscoelastic fluid. In a shear thinning fluid, we
find that the rotational speed is the fastest for prolate particles when they
are aligned with the shear-gradient direction, and slowest when they are
aligned with the flow direction. The converse is true for prolate particles in
a Boger fluid (e.g., second order fluid [16, 20]). This tendency is attributed
to the fact the hydrodynamic torque is due to shear stresses for the case of
shear thinning flows while the torque is due to normal stresses for the case
of viscoelastic flows. The different origins of the hydrodynamic torque thus
implies that the torque’s moment arm is the largest when the prolate particle
is oriented along the flow direction for viscoelastic fluids, while the moment
arm is the largest when the prolate particle is oriented along the shear
direction for shear thinning fluids. Similar trends are also seen for oblate
particles as well.

  * •

In a simple shear flow, the time period of tumbling of spheroids follow a non
monotonic trend with respect to Carreau number. At small Carreau number, the
time period increases from its Newtonian value, and then at high Carreau
number it decreases back to the same Newtonian value. On the other hand, in
the case of pressure driven flow, the time period of tumbling decreases
monotonically with Carreau number. Moreover, the time period in the different
Newtonian limits $(Cu\to 0)$ and $(Cu\to\infty)$ are different for pressure
driven flows. The time period in the Newtonian limit of high Carreau
$Cu\to\infty$ is $\beta$ times the time period in the Newtonian limit of low
Carreau ($Cu\to 0$) for pressure driven flows.

In the future, the present analysis may be extended to incorporate the effects
of other inelastic non Newtonian fluids like those with a spatially varying
viscosity [23, 50] or with spatially varying density [51, 52, 53]. Similarly,
instead of spheroids and ellipsoids, particles of other shapes like sheets
[54] , slender bodies [6] and non fore-aft symmetric particles [55] may be
analysed with the same theoretical framework. Other types of flow problems
beyond pressure driven flows, like sedimentation [18] for shear thinning
fluids may be analysed. Finally, higher order effects like those of inertia
[8, 6, 7], particle size [56], and effect of walls may also be taken into
account.

## Acknowledgments

The authors acknowledge funding from the American Chemical Society Petroleum
Research Fund (Grant No. ACS PRF 61266-DNI9).

## Appendix A: Derivation of reciprocal theorem for rigid body motion around
particle

Suppose we have two velocity fields in the same control volume $V$ of a
liquid. Both velocity fields satisfy the Stokes equations with a body force –
i.e.,

$\frac{\partial u_{i}}{\partial
x_{i}}=0;\qquad\frac{\partial\tau_{ij}^{N}}{\partial x_{j}}+b_{i}=0$ (35)

where $\tau_{ij}^{N}=\mu\dot{\gamma}_{ij}–p\delta_{ij}$ is the standard
Newtonian stress tensor with viscosity $\mu$, and $b_{i}$ is a spatially
varying body force. We will add superscripts to the symbols above to demarcate
the two flow fields – $(\alpha)$ is flow field one, while $(\chi)$ is flow
field two. These flow fields are related to each other via Green’s second
identity, which states that:

$\int_{S}u_{i}^{(\alpha)}\tau_{ij}^{N,(\chi)}n_{j}dS+\int_{V}u_{i}^{(\alpha)}b_{i}^{(\chi)}dV=\int_{S}u_{i}^{(\chi)}\tau_{ij}^{N,(\alpha)}n_{j}dS+\int_{V}u_{i}^{(\chi)}b_{i}^{(\alpha)}dV$
(36)

In the above expression, $S$ is the surface of the control volume and $n_{i}$
is the outward pointing normal vector for the control volume. While this
expression appears esoteric, it is quite powerful. It states that if one knows
information about one flow field (e.g., flow ($\chi$)), one can obtain
information about the other flow field (e.g., flow ($\alpha$)). This
expression can also be extended to non-Newtonian fluids, as will be
illustrated below.

We will let the two flow fields for our problem be the following:

##### Flow 1: Non-Newtonian flow around rigid particle

We will choose flow one (with superscript $\alpha$) to be the disturbance flow
around a particle with background velocity $u_{i}^{\infty}$ and body force
$b_{i}=\frac{\partial}{\partial x_{j}}\left(\tau_{ij}^{ex}\right)$. Thus, we
let:

$\displaystyle u_{i}^{(\alpha)}$ $\displaystyle=u_{i}–u_{i}^{\infty}$ (37a)
$\displaystyle\tau_{ij}^{N,(\alpha)}$
$\displaystyle=\tau_{ij}^{N}-\tau_{ij}^{N,\infty}$ (37b) $\displaystyle
b_{i}^{(\alpha)}$ $\displaystyle=\frac{\partial}{\partial
x_{j}}\left(\tau_{ij}^{ex}-\tau_{ij}^{ex,\infty}\right)$ (37c)

In the above expression, quantities with the superscript “$\infty$” denote
velocity and stress fields in the absence of the particle (i.e., due to
$u_{i}^{\infty}$ only), while quantities without the superscript are fields
with the particle present. We choose disturbance quantities for convenience,
as it will simplify the algebra later.

##### Flow 2: Stokes flow around particle from rigid body motion

We will let flow two (with superscript $\chi$) be Stokes flow around the same
particle undergoing rigid body motion. In other words, we let

$\displaystyle u_{i}^{(\chi)}$ $\displaystyle=v_{i}$ (38a)
$\displaystyle\tau_{ij}^{N,(\chi)}$ $\displaystyle=\Sigma_{ij}$ (38b)
$\displaystyle b_{i}^{(\chi)}$ $\displaystyle=0$ (38c)

where $v_{i}$ and $\Sigma_{ij}$ are the velocity and stress fields around the
particle from rigid body motion. On the particle surface,
$v_{i}=V_{i}+\epsilon_{ijk}\omega_{j}x_{k}$, where $x_{k}$ is the position
vector from the particle’s center of mass, while $V_{i}$ and $\omega_{i}$ are
the translational and rotational speeds. The external force and torque on the
particle will be $F_{i}^{aux}$ and $T_{i}^{aux}$. These will be related to the
translational and rotational speed through known resistance relationships.

Let us now substitute the information about the two flows into the integral
expression (36) above. We will choose the control volume $V$ to be the volume
outside of the particle. Since we are dealing with disturbance quantities, we
do not have to integrate over surfaces at infinity. We obtain:

$\int_{S_{p}}\left(u_{i}–u_{i}^{\infty}\right)\Sigma_{ij}n_{j}dS=\int_{S_{p}}v_{i}\left(\tau_{ij}^{N}-\tau_{ij}^{N,\infty}\right)n_{j}dS+\int_{V}v_{i}\frac{\partial}{\partial
x_{j}}\left(\tau_{ij}^{ex}-\tau_{ij}^{ex,\infty}\right)dV$ (39)

where $S_{p}$ is particle surface and $n_{i}$ is the normal vector pointing
into the particle (this is the outward pointing vector for the control volume
$V$). We can simplify the above expression using integration by parts on the
last integral, noting that the total stress tensor is
$\tau_{ij}=\tau_{ij}^{N}+\tau_{ij}^{ex}$. This yields:

$\int_{S_{p}}\left(u_{i}–u_{i}^{\infty}\right)\Sigma_{ij}n_{j}dS=\int_{S_{p}}v_{i}\left(\tau_{ij}-\tau_{ij}^{\infty}\right)n_{j}dS-\int_{V}\frac{\partial
v_{i}}{\partial x_{j}}\left(\tau_{ij}^{ex}-\tau_{ij}^{ex,\infty}\right)dV$
(40)

The next step in the derivation is to note that the velocity fields are rigid
body motion on the particle surface. Thus, on the surface,
$u_{i}=U_{i}+\epsilon_{ijk}\Omega_{j}x_{k}$, while for the other flow,
$v_{i}=V_{i}+\epsilon_{ijk}\omega_{j}x_{k}$. Substituting these expressions
gives:

$U_{i}F_{i}^{aux}+\Omega_{i}T_{i}^{aux}=V_{i}F_{i}^{ext}+\omega_{i}T_{i}^{ext}+\int_{S_{p}}u_{i}^{\infty}\Sigma_{ij}n_{j}dS-\int_{V}\frac{\partial
v_{i}}{\partial x_{j}}\left(\tau_{ij}^{ex}-\tau_{ij}^{ex,\infty}\right)dV$
(41)

In the above expression, $F_{i}^{ext}$ and $T_{i}^{ext}$ are the external
force and torque on the particle in the non-Newtonian flow (i.e., flow
$u_{i}$), while $F_{i}^{aux}$ and $T_{i}^{aux}$ are the external force and
torque on the particle from the other flow (i.e., flow $v_{i}$).

In the last part of the derivation, we note that all quantities associated
with the flow field $v_{i}$ are linear in the rigid body translation and
rotation ($V_{i}$, $\omega_{i}$). In other words, we can write the flow field
$v_{i}$ and stress field $\Sigma_{ij}$ in terms of these rigid body motions:

$\displaystyle v_{i}$
$\displaystyle=v_{ik}^{trans}V_{k}+v_{ik}^{rot}\omega_{k}$ (42a)
$\displaystyle\Sigma_{ij}$
$\displaystyle=\Sigma_{ijk}^{trans}V_{k}+\Sigma_{ijk}^{rot}\omega_{k}$ (42b)

In the above expression, $v_{ik}^{trans}$ is the flow field in the “$i$”
direction due to unit translation in the “$k$” direction. Similar notation
follows for the other quantities. For the force and torque $F_{i}^{aux}$ and
$T_{i}^{aux}$, we write them as:

$\displaystyle F_{i}^{aux}=R^{FU}_{ij}V_{j}+R^{F\Omega}_{ij}\omega_{j}$ (43a)
$\displaystyle T_{i}^{aux}=R^{TU}_{ij}V_{j}+R^{T\Omega}_{ij}\omega_{j}$ (43b)

where $R_{ij}^{FU}$, $R_{ij}^{F\Omega}$, $R_{ij}^{TU}$ and $R_{ij}^{T\Omega}$
are the Stokes-flow resistance tensors for the particle. If we substitute the
above two expressions (42) and (43) into the integral expression (41) and
perform some algebra, we obtain the final expression for the translational and
rotational velocity of the particle in a non-Newtonian fluid.

$\begin{bmatrix}R^{FU}_{ij}&R^{F\Omega}_{ij}\\\
R^{TU}_{ij}&R^{T\Omega}_{ij}\end{bmatrix}\begin{bmatrix}U_{j}\\\
\Omega_{j}\end{bmatrix}=\begin{bmatrix}F_{i}^{eff}\\\
T_{i}^{eff}\end{bmatrix}$ (44)

The effective force and torque are given by equations (20)-(22) in the main
text.

## Appendix B: Velocity fields from unit translation and rotation

From Kim and Karilla [4], the velocity fields $v_{ik}^{trans}$ and
$v_{ik}^{rot}$ are given by the following expressions. In these formulas, no
summation is assumed for repeated indices unless explicitly stated.

$\displaystyle
v_{ik}^{trans}=\frac{1}{16\pi\mu}R_{kk}^{FU}\left[\delta_{ik}G_{0}-x_{k}\frac{\partial
G_{0}}{\partial x_{i}}+\frac{a_{k}^{2}}{2}\frac{\partial^{2}G_{1}}{\partial
x_{i}\partial x_{k}}\right]$ (45a) $\displaystyle
v_{ik}^{rot}=\frac{3}{32\pi\mu}R_{kk}^{T\Omega}\sum_{j=1}^{3}\sum_{m=1}^{3}\epsilon_{jkm}\frac{\partial}{\partial
x_{m}}\left[\delta_{ij}G_{1}-x_{j}\frac{\partial G_{1}}{\partial
x_{i}}+\frac{a_{j}^{2}}{4}\frac{\partial^{2}G_{2}}{\partial x_{i}\partial
x_{j}}\right]$ (45b)

In the above expressions, the expression for $G_{n}$ is:

$G_{n}(x,y,z)=\int_{\lambda}^{\infty}\left(\frac{x^{2}}{a^{2}+t}+\frac{y^{2}}{b^{2}+t}+\frac{z^{2}}{c^{2}+t}-1\right)^{n}\frac{dt}{\Delta(t)}$
(46)

with $\Delta(t)=\sqrt{(a^{2}+t)(b^{2}+t)(c^{2}+t)}$ and $\lambda(x,y,z)$ being
the positive root of

$\frac{x^{2}}{a^{2}+t}+\frac{y^{2}}{b^{2}+t}+\frac{z^{2}}{c^{2}+t}=1$ (47)

## References

  * Purcell [1977] E. M. Purcell,  Life at low Reynolds number,  American Journal of Physics 45 (1977) 3–11. doi:10.1119/1.10903.
  * Happel and Brenner [1983] J. Happel, H. Brenner, Low Reynolds number hydrodynamics, Mechanics of fluids and transport processes, 2nd ed., Springer Netherlands, Dordrecht, 1983. doi:10.1007/978-94-009-8352-6.
  * Leal [2007] L. G. Leal, Advanced transport phenomena: fluid mechanics and convective transport processes, volume 7, Cambridge University Press, 2007\.
  * Kim and Karilla [2005] S. Kim, S. Karilla, Microhydrodynamics: Principles and Selected Application, 2 ed., Dover Publications, Mineola, New York, 2005.
  * Leal [1980] L. G. Leal,  PARTICLE MOTIONS IN A VISCOUS FLUID,  Ann. Rev. Fluid Mech 12 (1980) 435–76. URL: www.annualreviews.org.
  * Subramanian and Koch [2006] G. Subramanian, D. L. Koch,  Inertial effects on the orientation of nearly spherical particles in simple shear flow,  Journal of Fluid Mechanics 557 (2006) 257–296. doi:10.1017/S0022112006009724.
  * Dabade et al. [2016] V. Dabade, N. K. Marath, G. Subramanian,  The effect of inertia on the orientation dynamics of anisotropic particles in simple shear flow,  Journal of Fluid Mechanics 791 (2016) 631–703. doi:10.1017/jfm.2016.14.
  * Marath and Subramanian [2018] N. K. Marath, G. Subramanian,  The inertial orientation dynamics of anisotropic particles in planar linear flows,  Journal of Fluid Mechanics 844 (2018) 357–402. doi:10.1017/jfm.2018.184.
  * Einarsson et al. [2015a] J. Einarsson, F. Candelier, F. Lundell, J. R. Angilella, B. Mehlig,  Effect of weak fluid inertia upon Jeffery orbits,  Physical Review E - Statistical, Nonlinear, and Soft Matter Physics 91 (2015a). doi:10.1103/PhysRevE.91.041002.
  * Einarsson et al. [2015b] J. Einarsson, F. Candelier, F. Lundell, J. R. Angilella, B. Mehlig,  Rotation of a spheroid in a simple shear at small Reynolds number,  Physics of Fluids 27 (2015b). doi:10.1063/1.4921543.
  * D’Avino and Maffettone [2015] G. D’Avino, P. L. Maffettone,  Particle dynamics in viscoelastic liquids,  Journal of Non-Newtonian Fluid Mechanics 215 (2015) 80–104. doi:10.1016/j.jnnfm.2014.09.014.
  * D’Avino et al. [2017] G. D’Avino, F. Greco, P. L. Maffettone,  Particle Migration due to Viscoelasticity of the Suspending Liquid and Its Relevance in Microfluidic Devices,  Annual Review of Fluid Mechanics 49 (2017) 341–360. URL: https://doi.org/10.1146/annurev-fluid-010816-060150. doi:10.1146/annurev-fluid-010816-060150.
  * Einarsson et al. [2018] J. Einarsson, M. Yang, E. S. Shaqfeh,  Einstein viscosity with fluid elasticity,  Physical Review Fluids 3 (2018). doi:10.1103/PhysRevFluids.3.013301.
  * Murch and Shaqfeh [2020] W. L. Murch, E. S. Shaqfeh,  Collective effects in the sedimentation of particles in a viscoelastic fluid,  Physical Review Fluids 5 (2020). doi:10.1103/PhysRevFluids.5.073301.
  * Zhang et al. [2020] A. Zhang, W. L. Murch, J. Einarsson, E. S. G. Shaqfeh,  Lift and drag force on a spherical particle in a viscoelastic shear flow,  Journal of non Newtonian Fluid Mechanica (2020). URL: https://www.sciencedirect.com/science/article/pii/S0377025720300471.
  * Wang et al. [2020] S. Wang, C. W. Tai, V. Narsimhan,  Dynamics of spheroids in an unbound quadratic flow of a general second-order fluid,  Physics of Fluids 32 (2020). doi:10.1063/5.0030517.
  * Brunn [1980] P. Brunn,  THE MOTION OF RIGID PARTICLES IN VISCOELASTIC FLUIDS *,  Journal of Non-Newtonian Fluid Mechanicr 7 (1980) 271–288.
  * Kim [1986] S. Kim,  THE MOTION OF ELLIPSOIDS IN A SECOND ORDER FLUID,  Journal of Non-Newtonian Fluid Mechanics 21 (1986) 255–269.
  * Brunn [1977] P. Brunn,  The slow motion of a rigid particle in a second-order fluid,  Journal of Fluid Mechanics 82 (1977) 529–547. doi:10.1017/S0022112077000822.
  * Tai et al. [2020] C. W. Tai, S. Wang, V. Narsimhan,  Cross-stream migration of non-spherical particles in a second-order fluid-theories of particle dynamics in arbitrary quadratic flows,  Journal of Fluid Mechanics 895 (2020). doi:10.1017/jfm.2020.300.
  * Gong et al. [2022] J. Gong, V. A. Shaik, G. J. Elfring,  Active particles crossing sharp viscosity gradients (2022). URL: http://arxiv.org/abs/2208.09548. doi:https://doi.org/10.48550/arXiv.2208.09548.
  * Shaik and Elfring [2021] V. A. Shaik, G. J. Elfring,  On the hydrodynamics of active particles in viscosity gradients,  Physical Review of Fluids (2021).
  * Datt and Elfring [2019] C. Datt, G. J. Elfring,  Active Particles in Viscosity Gradients,  Physical Review Letters 123 (2019). doi:10.1103/PhysRevLett.123.158006.
  * Datt and Elfring [2018] C. Datt, G. J. Elfring,  Dynamics and rheology of particles in shear-thinning fluids,  Journal of Non-Newtonian Fluid Mechanics 262 (2018) 107–114. doi:10.1016/j.jnnfm.2018.03.016.
  * Abtahi and Elfring [2019] S. A. Abtahi, G. J. Elfring,  Jeffery orbits in shear-thinning fluids,  Physics of Fluids 31 (2019). doi:10.1063/1.5125468.
  * Pietrzyk et al. [2019] K. Pietrzyk, H. Nganguia, C. Datt, L. Zhu, G. J. Elfring, O. S. Pak,  Flow around a squirmer in a shear-thinning fluid,  Journal of Non-Newtonian Fluid Mechanics 268 (2019) 101–110. doi:10.1016/j.jnnfm.2019.04.005.
  * Elfring [2017] G. J. Elfring,  Force moments of an active particle in a complex fluid,  Journal of Fluid Mechanics 829 (2017). doi:10.1017/jfm.2017.632.
  * Van Gogh et al. [2022] B. Van Gogh, E. Demir, D. Palaniappan, O. S. Pak,  The effect of particle geometry on squirming through a shear-thinning fluid,  Journal of Fluid Mechanics 938 (2022). doi:10.1017/jfm.2022.116.
  * Elfring and Lauga [2015] G. J. Elfring, E. Lauga,  Theory of Locomotion Through Complex Fluids,  in: S. Spagnolie (Ed.), Complex Fluids in Biological Systems, Biological and Medical Physics, Biomedical Engineering, Springer Science+Business Media, 2015\. doi:10.1007/978-1-4939-2065-5{\\_}{\\_}8.
  * Datt et al. [2015] C. Datt, L. Zhu, G. J. Elfring, O. S. Pak,  Squirming through shear-thinning fluids,  Journal of Fluid Mechanics 784 (2015) R1. doi:10.1017/jfm.2015.600.
  * Nganguia et al. [2017] H. Nganguia, K. Pietrzyk, O. S. Pak,  Swimming efficiency in a shear-thinning fluid (2017). URL: http://arxiv.org/abs/1712.01305http://dx.doi.org/10.1103/PhysRevE.96.062606. doi:10.1103/PhysRevE.96.062606.
  * Qin and Pak [2023] K. Qin, O. S. Pak,  Purcell’s swimmer in a shear-thinning fluid,  Phys. Rev. Fluids 8 (2023) 33301. URL: https://link.aps.org/doi/10.1103/PhysRevFluids.8.033301. doi:10.1103/PhysRevFluids.8.033301.
  * Anand and Christov [2021] V. Anand, I. C. Christov,  Revisiting steady viscous flow of a generalized Newtonian fluid through a slender elastic tube using shell theory,  Z. Angew. Math. Mech. (ZAMM) 101 (2021) e201900309. doi:10.1002/zamm.201900309.
  * Anand and Christov [2019] V. Anand, I. C. Christov,  On the Deformation of a Hyperelastic Tube Due to Steady Viscous Flow Within,  in: H. Altenbach, A. Belyaev, V. A. Eremeyev, A. Krivtsov, A. V. Porubov (Eds.), Dynamical Processes in Generalized Continua and Structures, volume 103 of Springer Series on Advanced Structured Materials, Springer Nature, Cham, Switzerland, 2019, pp. 17–35. doi:10.1007/978-3-030-11665-1{\\_}2.
  * Anand et al. [2020] V. Anand, S. C. Muchandimath, I. C. Christov,  Hydrodynamic Bulge Testing: Materials Characterization Without Measuring Deformation,  ASME J. Appl. Mech. 87 (2020) 051012. doi:10.1115/1.4046297.
  * Anand et al. [2019] V. Anand, J. D. J. Rathinaraj, I. C. Christov,  Non-Newtonian fluid–structure interactions: Static response of a microchannel due to internal flow of a power-law fluid,  J. Non-Newtonian Fluid Mech. 264 (2019) 62–72. doi:10.1016/j.jnnfm.2018.12.008.
  * Venkatesh et al. [2022] A. Venkatesh, V. Anand, V. Narsimhan,  Peeling of linearly elastic sheets using complex fluids at low Reynolds numbers,  Journal of non Newtonian fluid mechanics (2022).
  * Anand and Christov [2020] V. Anand, I. C. Christov,  Transient compressible flow in a compliant viscoelastic tube,  Phys. Fluids 32 (2020) 112014. doi:10.1063/5.0022406.
  * Bird et al. [1987] R. B. Bird, R. C. Armstrong, O. Hassager, Dynamics of Polymeric Liquids, volume 1, 2 ed., John Wiley, New York, 1987\.
  * Anand [2016] V. Anand,  Effect of slip on heat transfer and entropy generation characteristics of simplified Phan-Thien-Tanner fluids with viscous dissipation under uniform heat flux boundary conditions: Exponential formulation,  Applied Thermal Engineering 98 (2016). doi:10.1016/j.applthermaleng.2015.12.025.
  * Anand [2015] V. Anand,  Entropy generation analysis of laminar flow of a nanofluid in a circular tube immersed in an isothermal external fluid,  Energy 93 (2015). doi:10.1016/j.energy.2015.09.019.
  * Chhabra and Richardson [2008] R. P. Chhabra, J. F. Richardson,  Non-Newtonian Fluid Behaviour,  in: Non-Newtonian Flow and Applied Rheology, 2 ed., Butterworth-Heinemann, Oxford, 2008, pp. 1–55. doi:10.1016/B978-0-7506-8532-0.00001-9.
  * Anand and Christov [2019] V. Anand, I. C. Christov,  On the enhancement of heat transfer and reduction of entropy generation by asymmetric slip in pressure-driven non-Newtonian microflows,  Journal of Heat Transfer 141 (2019). doi:10.1115/1.4042157.
  * Anand [2014] V. Anand,  Slip law effects on heat transfer and entropy generation of pressure driven flow of a power law fluid in a microchannel under uniform heat flux boundary condition,  Energy 76 (2014) 716–732. URL: http://dx.doi.org/10.1016/j.energy.2014.08.070. doi:10.1016/j.energy.2014.08.070.
  * Boyko and Stone [2021] E. Boyko, H. A. Stone,  Flow rate-pressure drop relation for shear-thinning fluids in narrow channels: Approximate solutions and comparison with experiments,  Journal of Fluid Mechanics 923 (2021). doi:10.1017/jfm.2021.621.
  * Pipe et al. [2008] C. J. Pipe, T. S. Majmudar, G. H. McKinley,  High shear rate viscometry,  in: Rheologica Acta, volume 47, 2008, pp. 621–642. doi:10.1007/s00397-008-0268-1.
  * Martin et al. [2019] C. P. Martin, S. Wang, S. Kim,  Surface tractions on an ellipsoid in Stokes flow: Quadratic ambient fields,  Physics of Fluids 31 (2019). doi:10.1063/1.5054698.
  * Wang et al. [2019] S. Wang, C. P. Martin, S. Kim,  Improper integrals as a puzzle for creeping flow around an ellipsoid,  Physics of Fluids 31 (2019). doi:10.1063/1.5050510.
  * Graham [2018] M. D. Graham, Microhydrodynamics, Brownian Motion, and Complex Fluids, 1 ed., Cambridge University Press, 2018. URL: www.cambridge.org/mathematics.
  * Dandekar and Ardekani [2020] R. Dandekar, A. M. Ardekani,  Swimming sheet in a viscosity-stratified fluid,  Journal of Fluid Mechanics (2020). doi:10.1017/jfm.2020.352.
  * Varanasi et al. [2022] A. K. Varanasi, N. K. Marath, G. Subramanian,  The rotation of a sedimenting spheroidal particle in a linearly stratified fluid,  Journal of Fluid Mechanics 933 (2022) A17. URL: https://www.cambridge.org/core/product/identifier/S0022112021010168/type/journal_article. doi:10.1017/jfm.2021.1016.
  * Varanasi and Subramanian [2022] A. K. Varanasi, G. Subramanian,  Motion of a sphere in a viscous density stratified fluid,  Journal of Fluid Mechanics 949 (2022) A29. URL: https://www.cambridge.org/core/product/identifier/S0022112022007601/type/journal_article. doi:10.1017/jfm.2022.760.
  * More and Ardekani [2022] R. V. More, A. M. Ardekani,  Motion in Stratified Fluids,  Annual Review of Fluid Mechanics (2022). URL: https://doi.org/10.1146/annurev-fluid-120720-011132. doi:10.1146/annurev-fluid-120720-011132.
  * Dandekar et al. [2019] R. Dandekar, V. A. Shaik, A. M. Ardekani,  Swimming sheet in a density-stratified fluid,  Journal of Fluid Mechanics 874 (2019) 210–234. doi:10.1017/jfm.2019.445.
  * Candelier and Mehlig [2016] F. Candelier, B. Mehlig,  Settling of an asymmetric dumbbell in a quiescent fluid,  Journal of Fluid Mechanics 802 (2016) 174–185. doi:10.1017/jfm.2016.350.
  * Anand and Subramanian [2023] P. Anand, G. Subramanian,  Inertial migration in pressure-driven channel flow: beyond the Segre-Silberberg pinch (2023).