Electrokinetically enhanced cross-stream particle migration in viscoelastic flows

TL;DR

This paper explores how electrokinetic effects can significantly enhance particle migration in viscoelastic flows within microfluidic devices, enabling improved particle separation based on electrical properties.

Contribution

It introduces a combined experimental and theoretical study showing electrokinetics can amplify particle migration and control particle focusing in viscoelastic flows.

Findings

Electrokinetics can enhance particle migration by an order of magnitude.

Particles can be focused at different positions based on electrical properties.

Theoretical predictions agree with experimental results.

Abstract

Advancements in understanding the lateral migration of particles have helped in enhanced separation in microfluidic devices. In this work, we investigate the effects of electrokinetics on the particle migration in a viscoelastic flow, where the electric field is applied parallel to the flow. Through experiments and theoretical analysis, we show that the interaction of electrokinetic and rheological effects can result in an enhancement in the migration by an order of magnitude. The theoretical analysis, in agreement with experiments, demonstrates that the particles can be focused at different equilibrium positions based on their intrinsic electrical properties.

Figures

Equations

$${\nabla }^2\varphi =0$$

$${\mathbf{e}}_r\text{bcdot}\text{bnabla}\varphi =0\text{mbox}atr=1\text{mbox}and{\mathbf{e}}_z\text{bcdot}\text{bnabla}\varphi =0\text{mbox}atz=-s/\kappa \&z=(1-s)/\kappa .$$

$$\varphi \to {\varphi }_{\infty }\text{mbox}asr\to \infty .$$

$$\text{bnabla}\text{bcdot}\mathbf{U}=0\text{mbox}and{\text{Rey}}_p\left(\mathbf{U}\text{bcdot}\text{bnabla}\mathbf{U}\right)=\text{bnabla}\text{bcdot}\text{mathsfbi}\mathbf{T},$$

$$\text{mathsfbi}\mathbf{T}=-P\text{mathsfbi}\mathbf{I}+{\text{mathsfbi}\mathbf{E}}^{(1)}+De\text{mathsfbi}\mathbf{S}.$$

$$\text{mathsfbi}\mathbf{S}={\text{mathsfbi}\mathbf{E}}^{(1)}\text{bcdot}{\text{mathsfbi}\mathbf{E}}^{(1)}+\delta {\text{mathsfbi}\mathbf{E}}^{(2)}.$$

$${\text{mathsfbi}\mathbf{E}}^{(2)}=\mathbf{U}\text{bcdot}\text{bnabla}{\text{mathsfbi}\mathbf{E}}^{(1)}+{\text{mathsfbi}\mathbf{E}}^{(1)}\text{bcdot}\text{bnabla}{\mathbf{U}}^{†}+\text{bnabla}\mathbf{U}\text{bcdot}{\text{mathsfbi}\mathbf{E}}^{(1)},$$

$$De=\frac{({\Psi }_1^{\ast }+{\Psi }_2^{\ast })U_{max}^{\ast }}{{\mu }^{\ast }{l}^{\ast }},\delta =\frac{-{\Psi }_1^{\ast }}{2({\Psi }_1^{\ast }+{\Psi }_2^{\ast })},{\Psi }_1^{\ast }=\frac{T_{xx}^{\ast }-T_{zz}^{\ast }}{{\left(\partial U_x^{\ast }/\partial z^{\ast }\right)}^2},{\Psi }_2^{\ast }=\frac{T_{yy}^{\ast }-T_{zz}^{\ast }}{{\left(\partial U_x^{\ast }/\partial z^{\ast }\right)}^2}.$$

$$\mathbf{U}={\mathbf{\Omega }}_s\times \mathbf{r}+Ha{\zeta }_p\nabla \varphi \text{mbox}atr=1.$$

$$\mathbf{U}=Ha{\zeta }_w\text{bnabla}\varphi -{\mathbf{U}}_s\text{mbox}atwalls.$$

$$\text{mathsfbi}\mathbf{M}=\text{bnabla}\varphi \text{bnabla}\varphi -\frac{1}{2}\left(\text{bnabla}\varphi \text{bcdot}\text{bnabla}\varphi \right)\text{mathsfbi}\mathbf{I}.$$

$$\mathbf{F}={\mathbf{F}}_H+{\mathbf{F}}_M=\underset{S_p}{\int }\mathbf{n}\text{bcdot}\text{mathsfbi}\mathbf{T}\mathrm{d}S+4\pi Ha\underset{S_p}{\int }\mathbf{n}\text{bcdot}\text{mathsfbi}\mathbf{M}\mathrm{d}S=0,$$

$$\mathbf{L}={\mathbf{L}}_H+{\mathbf{L}}_M=\underset{S_p}{\int }\mathbf{n}\times \left(\mathbf{n}\text{bcdot}\text{mathsfbi}\mathbf{T}\right)\mathrm{d}S+4\pi Ha\underset{S_p}{\int }\mathbf{n}\times \left(\mathbf{n}\text{bcdot}\text{mathsfbi}\mathbf{M}\right)\mathrm{d}S=0.$$

$$\mathbf{U}{}_{\infty }=\left(\alpha +\beta z+\gamma z^2\right){\mathbf{e}}_x+Ha{\zeta }_w\text{bnabla}{\varphi }_{\infty }-{\mathbf{U}}_s.$$

$$P{}_{\infty }=2\gamma x+4De(1+2\delta )(\beta z+\gamma z^2)+\mathcal{I},$$

$$\alpha =4s\left(1-s\right)/\kappa ,\beta =4\left(1-2s\right),\gamma =-4\kappa ,$$

$${\nabla }^2\psi =0.$$

$${\mathbf{e}}_{\mathbf{r}}\text{bcdot}\text{bnabla}\psi =-{\mathbf{e}}_{\mathbf{r}}\text{bcdot}\text{bnabla}{\varphi }_{\infty }\text{mbox}atr=1\&{\mathbf{e}}_{\mathbf{z}}\text{bcdot}\text{bnabla}\psi =0\text{mbox}atwalls.$$

$$\text{bnabla}\text{bcdot}\mathbf{v}=0\text{mbox}and\text{bnabla}\text{bcdot}{\text{mathsfbi}\mathbf{\sigma }}_H=0.$$

$$\text{mathsfbi}\mathbf{s}={\text{mathsfbi}\mathbf{e}}^{(1)}\text{bcdot}{\text{mathsfbi}\mathbf{e}}^{(1)}+{\text{mathsfbi}\mathbf{w}}^{(1)}+\delta ({\text{mathsfbi}\mathbf{e}}^{(2)}+{\text{mathsfbi}\mathbf{w}}^{(2)}),$$

$${\text{mathsfbi}\mathbf{w}}^{(1)}$$

$${\text{mathsfbi}\mathbf{w}}^{(2)}$$

$$\mathbf{v}={\mathbf{\Omega }}_s\times \mathbf{r}+Ha{\zeta }_p\text{bnabla}\left(\psi +{\varphi }_{\infty }\right)-\mathbf{U}{}_{\infty }\text{mbox}atr=1,$$

$$\mathbf{v}=Ha{\zeta }_w\text{bnabla}\psi \text{mbox}atwalls,$$

$$\mathbf{v}\to 0\text{mbox}asr\to \infty .$$

$$\text{mathsfbi}\mathbf{M}$$

$${\text{mathsfbi}\mathbf{\sigma }}_M=\text{bnabla}\psi \text{bnabla}\psi -\frac{1}{2}\left(\text{bnabla}\psi \text{bcdot}\text{bnabla}\psi \right)\text{mathsfbi}\mathbf{I}.$$

$$\mathbf{F}={\mathbf{F}}_H+{\mathbf{F}}_M=\underset{S_p}{\int }\mathbf{n}\text{bcdot}{\text{mathsfbi}\mathbf{\sigma }}_H\mathrm{d}S+4\pi Ha\underset{S_p}{\int }\mathbf{n}\text{bcdot}{\text{mathsfbi}\mathbf{\sigma }}_M\mathrm{d}S=0,$$

$$\mathbf{L}={\mathbf{L}}_H+{\mathbf{L}}_M=\underset{S_p}{\int }\mathbf{n}\times \left(\mathbf{n}\text{bcdot}{\text{mathsfbi}\mathbf{\sigma }}_H\right)\mathrm{d}S+4\pi Ha\underset{S_p}{\int }\mathbf{n}\times \left(\mathbf{n}\text{bcdot}{\text{mathsfbi}\mathbf{\sigma }}_M\right)\mathrm{d}S=0.$$

$${\mathbf{F}}_M=4\pi Ha\left(\frac{3\pi }{16}{\kappa }^4\left(\mathfrak{Z}\left(4,s\right)-\mathfrak{Z}\left(4,1-s\right)\right)\right){\mathbf{e}}_z+O\left({\kappa }^7\right)\text{and}{\mathbf{L}}_M=0.$$

$$U_{mig}^M\approx \frac{4\pi Ha}{6\pi (1+O(\kappa ))}\left(\frac{3\pi }{16}{\kappa }^4\left(\mathfrak{Z}\left(4,s\right)-\mathfrak{Z}\left(4,1-s\right)\right)\right).$$

$$\xi ={\xi }_{(0)}+De{\xi }_{(1)}+\cdots .$$

$$\left.\begin{array}[]{l}\qquad\;\quad\bnabla\bcdot{\boldsymbol{v}_{(0)}}=0,\\ {\nabla^{2}}{\boldsymbol{v}_{(0)}}-\bnabla{p_{(0)}}=\boldsymbol{0},\\ \qquad\qquad\quad{\boldsymbol{v}_{(0)}}={\boldsymbol{\Omega}_{s}}_{(0)}\times\boldsymbol{r}+H\!a{\zeta_{p}}\bnabla\left({\psi+{\phi_{\infty}}}\right)-{\boldsymbol{U}\!\!_{\infty}}_{(0)}\quad\mbox{at\ }r=1,\\ \qquad\qquad\quad{\boldsymbol{v}_{(0)}}=H\!a{\zeta_{w}}\bnabla\psi\quad\mbox{at walls},\\ \qquad\qquad\quad{\boldsymbol{v}_{(0)}}\to\boldsymbol{0}\quad\mbox{as\ }r\to\infty,\end{array}\right\}$$

Full text

1 Introduction

In the past two decades, biological and healthcare research has significantly benefited by the application of various particle and cell sorting techniques—inertial focusing (Di Carlo, 2009; Warkiani et al., 2015; Zhang et al., 2016) and viscoelastic focusing (Leshansky et al., 2007; Lu et al., 2017)—in microfluidic devices. The working principle of such devices is based on the exploitation of non-linearity in the flow (arising from inertial or polymeric stresses) which breaks the symmetry and generates a lateral force on the particle, pushing it towards the equilibrium positions. Since most biological fluids are non-Newtonian, it is essential to understand the rheological effects on particle migration. The pioneering studies of Mason and his coworkers (Karnis & Mason, 1966; Gauthier et al., 1971) and Ho & Leal (1976) clearly demonstrated that the presence of polymer chains in the fluid exert a normal stress on the particle, pushing it towards the region of lowest shear. In biological fluids, strands of DNA enact the role of polymeric chains and help in the lateral migration of cells (Kim & Kim, 2016). In recent years, this normal stress induced migration has been exploited to design tubular-viscoelastic-focusers, in which the region of lowest shear is the channel axis (D’Avino et al., 2017). However, the labor-intensive fabrication process limits the accessibility of these devices (Nam et al., 2015). The more readily available rectangular microchannels are plagued by the limitation of multiple equilibrium positions (1 center and 4 corners), which renders the downstream detection and separation difficult (Lu et al., 2017; Tian et al., 2019). Furthermore, the passive nature of the migration requires large focusing lengths which are dependent on the flowrate (D’Avino et al., 2017).

Unlike flows at large scale, in micro-scale flows, the short-ranged effects near the particle surface can determine their fate (Datt et al., 2017). One can resort to manipulation of surface effects through external fields to tune particle focusing (Rossi et al., 2019). This has potential applications in focusing, separation and detection of rare circulating tumor cells in the blood stream (Zhang et al., 2014). Recent studies have shown that cancer cells exhibit highly negative charge on their surface which can be used to separate them from the blood stream by applying external electric fields (Chen et al., 2016; Le et al., 2019). In our previous study (Choudhary et al., 2019), we studied the effects of electric fields on the inertial migration in pressure driven flow of a Newtonian fluid. Since most biological fluids exhibit non-Newtonian behaviour, it is important to understand the influence of electrokinetics on the migration in complex fluids.

In this work, we experimentally and theoretically investigate the influence of electrophoresis on viscoelastic focusing of particles, where the electric field is applied parallel to the flow. A vast majority of cells and particles, when submerged in aqueous electrolytes, acquire a surface charge which attracts the counter-charged ions. As a result a thin double layer forms around the particle (Anderson, 1989). An external electric field causes a tangential slip on the surface of the particle. This generates a disturbance field around the particle which may amplify the asymmetry in normal stresses. Our experimental results demonstrate that the focusing is significantly accelerated in the presence of an electric field, which reduces the overall focusing length. Furthermore, the electric field provides an external control over the particle migration and enables charge-selective focusing. We explain the experimental observations qualitatively and quantitatively using an approach based on perturbation theory: we analyze the system in the limit of slow and slowly varying flows and capture the viscoelastic effects using a second order fluid model (Ho & Leal, 1976; Bird et al., 1987). This allows us to predict the migration velocity through a perturbation solution. Using the reciprocal theorem we show that the electrophoretic motion of the particle generates additional polymeric stresses which enhance the migration by $O(\kappa^{-1})$, where $\kappa$ is the particle to channel size ratio. The theoretical results (obtained without having to resort to fit any parameter) are found to be in reasonable agreement with the experimental observations.

The rest of this article is organized as follows. In §2 we describe the experimental conditions and results. The explanation of these results is conducted through a mathematical model which is elaborated in §3. In §4 we qualitatively and quantitatively compare the theoretical and experimental results. Key conclusions and potential future directions are described in §5.

2 Experiments

In our experiments, we use a 2 cm long PDMS microchannel with a rectangular cross-section of $66\times 54$ $\mu$m2. The electric fields were generated using a high-voltage DC power supply using a wire connection to electrodes at the channel inlet and outlet. We dissolved polyethylene oxide powder (PEO, molecular weight = 2$\times 10^{6}$ g/mol, Sigma-Aldrich) into a buffer solution (0.01 mM phosphate buffer mixed with 0.5$\%$ Tween 20, Fisher Scientific) at a concentration of 250 ppm. A dilute suspension ($<$0.1$\%$ in volume fraction) of polystyrene spheres of 2.2 $\mu$m diameter is prepared in the PEO solution and is pumped at 25 $\rm{\mu L/h}$. The particle zeta potential was measured to be $-83$ mV. Particle motion was visualized at the outlet of the microchannel through an inverted microscope (Nikon Eclipse TE2000U) equipped with a CCD camera (Nikon DS-Qi1MC). Digital videos were recorded at a rate of 15 frames per second, from which the superimposed images were obtained and further processed in the Nikon imaging software (NIS-Elements AR 3.22).

Fig.1(a) shows the particles at the inlet and outlet of the channel when no electric field is applied i.e., the focusing is purely viscoelastic (Karnis & Mason, 1966; Gauthier et al., 1971; Ho & Leal, 1976). The particles remain almost uniformly distributed at the channel exit. This weak focusing is depicted in Fig.1(d), which shows the change in particle stream width (lateral width of the particle distribution) along the channel length. In Fig.1(b) an application of electric field parallel to the flow shows that a significant fraction of the particles are focused at the centerline. Fig.1(d) shows that the introduction of electric field results in substantial reduction in particle stream width. Upon reversing the direction of the electric field (see Fig.1(c)), the particles focus near the walls. These results demonstrate the potential to efficiently separate particles and cells based on their intrinsic electrical properties. To extrapolate these observations towards potential applications, it is pivotal to understand the role of electrical properties and flow profile characteristics.

In what follows, we formulate a mathematical framework to predict, quantify and draw physical insights into the migration behavior observed in the experiments. Since the particle suspension in the experimental study is dilute, we analyze the migration of a single particle. We model the fluid medium using a second-order fluid (SOF) model (Bird et al., 1987), as we are only interested in capturing the effect of small deviations from the Newtonian behavior. Although the SOF model is only valid for slow and slowly varying (thus nearly Newtonian) flows, it’s an asymptotic approximation for most viscoelastic flows (Ho & Leal, 1976). Furthermore, it has been shown in the past that dilute PEO solutions ($C_{PEO}^{*}<500$ ppm) exhibit pure elastic behavior (i.e. negligible shear thinning) (Rodd et al., 2005; Yang et al., 2011; Kim et al., 2012). Hence employing the SOF model can help us capture the dynamics in slow Boger fluid flows, both qualitatively and quantitatively.

\floatsetup

[figure]style=plain

3 Mathematical model

Figure 2 shows the schematic of a neutrally buoyant electrophoretic sphere suspended in a pressure driven flow of SOF, subjected to a parallel electric field ($\mathcal{E}_{\infty}^{*}$). The fluid moves with a maximum velocity of $U_{max}^{*}$ at the centerline. The frame of reference is placed at the center of the particle. Typical PDMS micro-channel walls and polystyrene particles usually acquire a surface charge (formally quantified as zeta-potential $\zeta^{*}$) when in contact with an electrolytic solution (Sze et al., 2003), which gives rise to an electrical double layer (EDL) (Anderson, 1989). In comparison to the channel and particle size, the EDL can be assumed to be asymptotically thin i.e., $\lambda^{*}_{D}/a^{*}\ll 1$ (where $\lambda^{*}_{D}$ is the EDL thickness). Under the influence of an external electric field, the ions inside the EDL slip in the tangential direction, generating an electro-osmotic slip. For low zeta potentials ($\zeta_{p}^{*}\sim 25$ mV), the slip at the particle surface can be represented by Smoluchowski’s relation (Keh & Anderson, 1985; Anderson, 1989). Recently, Ghosh et al. (2016) showed that the surface slip in the SOF does not get modified for a surface with constant zeta potential. Hence we employ the Smoluchowski slip in this work.

3.1 Governing equations

We represent the system in terms of dimensionless variables, where the coordinates, electrostatic potential, velocity and pressure are non-dimensionalized using $a^{*}$, $|\mathcal{E}_{\infty}^{*}|a^{*}$, $\kappa U_{max}^{*}$ and $\mu^{*}\kappa U_{max}^{*}/a^{*}$, respectively ($\mu^{*}$ being the dynamic viscosity). Here $\kappa$ is the ratio of particle radius ($a^{*}$) to the channel width ($l^{*}$). To obtain analytical results, we will later restrict our analysis to small particles ($\kappa\ll 1$), which is consistent with the experimental conditions.

The potential ($\phi$) in the electro-neutral bulk region is governed by the Laplace equation

[TABLE]

The boundary conditions at the particle and wall surface (i.e. slip plane of the EDL) are governed by the no-flux condition

[TABLE]

where $s\equiv d^{*}/l^{*}$ (as defined in fig. 2). The potential is undisturbed ($\phi_{\infty}\sim-x$) far away from the particle:

[TABLE]

The hydrodynamics is governed by the continuity and the Navier-Stokes equations:

[TABLE]

where $\Rey_{p}$ is particle Reynolds number defined as: ${\rho^{*}U_{\max}^{*}\kappa a^{*}}/{\mu^{*}}$ ($\rho^{*}$ being the fluid density). We model the hydrodynamic stress using the second-order-fluid model (Bird et al., 1987):

[TABLE]

Here, $\boldsymbol{\mathsfbi{E}}^{(1)}$ represents the rate of strain tensor ($\bnabla\boldsymbol{U}+\bnabla\boldsymbol{U}^{\dagger}$), $\dagger$ denotes the transpose, $\boldsymbol{\mathsfbi{S}}$ represents the polymeric stress:

[TABLE]

$\boldsymbol{\mathsfbi{E}}^{(2)}$ is the steady-state covariant convected derivative of $\boldsymbol{\mathsfbi{E}}^{(1)}$ (also known as the Rivlin-Ericksen tensor):

[TABLE]

Here, $De$ (Deborah number) is the ratio of viscoelastic time scale ($(\Psi_{1}^{*}+\Psi_{2}^{*})/\mu^{*}$) to the time scale corresponding to the average shear in the background flow ($l^{*}/U_{max}^{*}$), $\delta$ is a viscometric parameter, $\Psi_{1}^{*}$ and $\Psi_{2}^{*}$ are the steady shear viscometric coefficients (Bird et al., 1987):

[TABLE]

Since the migration time scale is much larger than the convective time scale ($a^{*}/\kappa U_{max}^{*}$), we analyze the migration in quasi steady-state (Becker et al., 1996; Choudhary et al., 2019). Also, we neglect inertial effects (i.e. $\Rey_{p}\ll De$) which leaves the polymeric stress as the sole non-linearity in the system.

The sphere moves with a translational ($\boldsymbol{U}_{s}$) and rotational ($\boldsymbol{\Omega}_{s}$) velocity, which is to be determined from force and torque balance. Since the domain is infinite in the y-direction, we can consider: ${\boldsymbol{U}_{s}}={U_{sx}}{\boldsymbol{e}_{x}}+{U_{sz}}{\boldsymbol{e}_{z}}$ and ${\boldsymbol{\Omega}_{s}}={\Omega_{sy}}{\boldsymbol{e}_{y}}$. To account for the electrokinetic effects in the EDL, we employ the slip-boundary condition derived by Smoluchowski (1903) and impose at the particle surface:

[TABLE]

Here, $\zeta_{p}$ is the dimensionless particle surface zeta potential ($\zeta_{p}^{*}/|\mathcal{E}_{\infty}^{*}|a^{*}$). $H\!a$ is a Hartmann-type number: $H\!a=\varepsilon^{*}{\mathcal{E}^{*}}_{\infty}^{2}a^{*}/4\pi\mu^{*}\kappa U_{max}^{*}$, which denotes the ratio of electrical energy density to shear stress. Here $\varepsilon^{*}$ is the electrical permittivity of the medium. Similarly, the flow satisfies the Smoluchowski slip condition at the walls

[TABLE]

As the particle is freely suspended, the translational ($\boldsymbol{U}_{s}$) and rotational velocity ($\boldsymbol{\Omega}_{s}$) of the particle are obtained using force-free and torque-free conditions: the total force and torque on the particle must be zero. Apart from the hydrodynamic stress, we also have to account for the Maxwell stress in the calculation of total force and torque. Maxwell stress arises due to a lateral asymmetry in the electrostatic potential distribution around the particle due to the presence of the wall. The Maxwell stress is defined as

[TABLE]

It should be noted that there is no contribution from Maxwell stress to the fluid momentum equation (4) because $\bnabla\bcdot\boldsymbol{\mathsfbi{M}}=0$ in the electroneutral fluid (Yariv, 2006). The total force and torque free conditions for particle are given by

[TABLE]

Here subscript $H$ and $M$ represent contributions from hydrodynamic and Maxwell stress, respectively.

3.1.1 Undisturbed flow

In the absence of particle, the flow profile for a fully developed pressure-driven second-order fluid flow (in the frame of reference of the particle moving with $\boldsymbol{U}_{s}$) is

[TABLE]

This is identical to that obtained for a Newtonian fluid (Choudhary et al., 2019). However, the pressure exhibits a variation in the lateral direction for SOF

[TABLE]

where $\mathcal{I}$ is an integration constant. The constants $\alpha,\beta$ and $\gamma$ are:

[TABLE]

where $\beta$ and $\gamma$ represent the shear and curvature of the background flow.

3.2 Equations governing the disturbance

We analyze the system in terms of disturbance variables (i.e. deviation from the undisturbed flow). For the different variables, these are defined as: potential ($\psi=\phi-\phi_{\infty}$), velocity ($\boldsymbol{v}=\boldsymbol{U}-\boldsymbol{U}\!\!_{\infty}$) and pressure ($p=P-P_{\infty}$). Since the bulk fluid medium is electroneutral, the disturbance potential is also governed by the Laplace equation

[TABLE]

Far away from the particle, $\psi$ decays to zero. Using the definition of disturbance potential in (2):

[TABLE]

The disturbance flow is governed by continuity and Cauchy’s momentum equation:

[TABLE]

Here $\boldsymbol{\sigma}_{H}$ is the hydrodynamic disturbance stress tensor which contains the Newtonian and polymeric contributions ($\boldsymbol{\mathsfbi{\sigma}}_{H}=-p\boldsymbol{\mathsfbi{I}}+\boldsymbol{\mathsfbi{e}}^{(1)}+De\,\boldsymbol{\mathsfbi{s}}$). Here, $\boldsymbol{\mathsfbi{I}}$ is the identity tensor, $\boldsymbol{\mathsfbi{e}}^{(1)}$ is the rate of strain tensor for the disturbance flow ($\bnabla\boldsymbol{v}+\bnabla\boldsymbol{v}^{\dagger}$), and $\boldsymbol{\mathsfbi{s}}$ is the disturbance polymeric stress tensor

[TABLE]

[TABLE]

Here, $\boldsymbol{\mathsfbi{E}}_{\infty}^{(1)}$ is the rate of strain tensor for the undisturbed flow; $\boldsymbol{\mathsfbi{e}}^{(2)}$ is the disturbance Rivlin-Eriksen (RE) tensor; $\boldsymbol{\mathsfbi{w}}^{(1)}$ is the ‘interaction tensor’ (arising from the interaction between background flow and disturbance flow); and $\boldsymbol{\mathsfbi{w}}^{(2)}$ is the RE interaction tensor.

Using (9) and (10), we find that the disturbance flow is subject to the following boundary conditions:

[TABLE]

The total Maxwell stress can be decomposed into four components as

[TABLE]

The first term arises from the undisturbed component, the second and third terms denote the contribution from the interaction between applied ($\phi_{\infty}$) and disturbed potential field ($\psi$), and the fourth term denotes the disturbance Maxwell stress. Yariv (2006) showed that the first three components do not result in a net force or torque and the leading order contribution arises from the fourth term111Since $\boldsymbol{\mathsfbi{M}}$ is divergence free in the fluid domain, only the terms (in eq.22) decaying as $O(r^{-2})$ contribute to the force (Guazzelli & Morris, 2011). Only the disturbance Maxwell stress generates $O(r^{-2})$ terms and provides the leading order contribution (Yariv, 2006, p.3). . Thus, the disturbance Maxwell stress is denoted as

[TABLE]

Similarly, the undisturbed hydrodynamic field cannot generate a net force or torque on the particle. Hence, the force and torque free conditions are given by

[TABLE]

3.3 Migration due to Maxwell stress

Here we evaluate the force and lateral migration arising from the Maxwell stress (23) alone. Since the potential distribution is not affected by the fluid rheology, the Maxwell force on the particle in a second order fluid is identical to that in a Newtonian fluid (Yariv, 2006). Following our previous analysis (Choudhary et al., 2019), we use the method of reflections to determine the potential distribution in the particle-wall configuration (Brenner, 1962). Since our aim is to obtain qualitative insights and analytical results, we restrict the analysis to small particles i.e. $\kappa\ll 1$. Assuming the small particle is not too close to the walls (i.e. $s\gg\kappa$), the disturbance potential is sought as successive reflections: $\psi={{}_{(1)}\psi}+{{}_{(2)}\psi}+\dotsm$. Here, ${}_{(i)}\psi$ represents the $i^{th}$ reflection, where the odd reflections satisfy boundary conditions at the particle surface and even reflections satisfy the wall boundary conditions.

Accounting for each successive pair of reflections increases the accuracy by $O(\kappa)$ (Happel & Brenner, 2012, p.286). Since $\kappa\ll 1$, we take only the first two reflections of disturbance potential into account. Evaluation of $\boldsymbol{F}_{M}$ and $\boldsymbol{L}_{M}$ (detailed in Appendix A.2) yields:

[TABLE]

The Maxwell force acts along the cross-stream direction alone and is wall-repulsive in nature; the Maxwell torque is zero. The contribution of Maxwell stresses to the migration velocity is

[TABLE]

Here, superscript $M$ shows. that the contribution is from Maxwell stress, $\mathfrak{Z}$ is the generalized Riemann zeta function, and $6\pi(1+O(\kappa))$ is the viscous resistance which includes wall correction at the first order (Brenner, 1961, p.246)222The details of evaluation through the method of reflections are provided in the supplementary material.. Later we show that $U_{mig}^{M}$ decays rapidly away from the walls and is negligible in the bulk of the channel, where migration due to hydrodynamic stress dominates.

3.4 Migration due to hydrodynamic stress

3.4.1 Perturbation expansion

We find the migration at $O(De)$ using a regular perturbation expansion. For small values of $De$, the disturbance field variables are expanded as:

[TABLE]

Here, $\xi$ is a generic field variable which represents velocity ($\boldsymbol{v}$), pressure ($p$), translational ($\boldsymbol{U}_{s}$) and angular velocity ($\boldsymbol{\Omega}_{s}$). Since the electrostatic potential is decoupled from the hydrodynamics, the variable $\psi$ is not expanded as in (27). We substitute (27) in the equations governing the hydrodynamics (18-21). We obtain the problem at $O(1)$ as

[TABLE]

and at $O(De)$ as:

[TABLE]

Here, $\boldsymbol{\mathsfbi{s}}_{(0)}$ is the disturbance polymeric stress of $O(De^{0})$ field. The angular ($\boldsymbol{\Omega}_{s}$) and translational velocity ($\boldsymbol{U_{s}}$) are evaluated by using (24), at each order.

3.4.2 Reciprocal theorem

The symmetry of Stokes flow implies that the $O(1)$ field cannot produce a lateral lift. Therefore, the lift must arise from $O(De)$ field. The reciprocal theorem (Ho & Leal, 1976) allows us to determine the migration associated with $\boldsymbol{v}_{(1)}$, without having to solve for $O(De)$ field. The reciprocal theorem relates the properties of an unknown Stokes flow to a known test flow field, provided both flow fields correspond to the same geometry. The test field ($\boldsymbol{u}^{t},\,p^{t}$) is chosen to be that generated by a moving sphere in the positive z-direction (towards the upper wall) with unit velocity in a quiescent Newtonian fluid (Refer Appendix B for the details of test field). Following (Ho & Leal, 1976, p.789), we obtain the lateral migration velocity at $O(De)$ as:

[TABLE]

Here, superscript $H$ denotes the hydrodynamic contribution and $6\pi(1+O(\kappa))$ is the viscous resistance including wall correction at the first order. The details of derivation are shown in the supplementary material.

3.4.3 Evaluation of $\boldsymbol{v}_{(0)}$

Following the method used for the evaluation of potential ($\psi$), we use the method of reflections to find the velocity field which satisfies the boundary conditions at the particle and wall surface with sufficient accuracy (Brenner, 1962). As in §3.3, we restrict the analysis to small particles i.e. $\kappa\ll 1$. This condition is in addition to $De\ll 1$, which was used to expand the field variables in (27). We seek the disturbance ($\boldsymbol{v}_{(0)},p_{(0)}$) as

[TABLE]

Here ${}_{(i)}\boldsymbol{v}_{\left(0\right)}$ represents the $i$th reflection of the creeping flow disturbance velocity. Since $\kappa\ll 1$, we take the first two reflections into account. We substitute the above expansion in the equations governing the creeping flow hydrodynamics (28) and obtain the following set of problems:

[TABLE]

[TABLE]

The first reflection (${}_{(1)}\boldsymbol{v}_{(0)},\,{}_{(1)}p_{(0)}$) is found by employing Lamb’s general solution (Lamb, 1975; Happel & Brenner, 2012):

[TABLE]

Here, the coefficients are defined as:

[TABLE]

The terms multiplying the coefficients $A_{1}$, $B_{1}$, $C_{1}$, $D_{1}$, $E_{1}$ represent the force-monopole (stokeslet), source-dipole, antisymmetric dipole (rotlet), symmetric dipole (stresslet), and octupole singularities, respectively (Kim & Karrila, 2013; Guazzelli & Morris, 2011).

The other disturbances (terms multiplying $F_{1},\,G_{1},\,H_{1}$) are further singularities in the multipole expansion which arise due to the curvature ($\gamma$) in the background flow field.

The evaluation of higher order reflections can be found using Faxén transformation (Faxén, 1922) and is identical to our previous analysis (Choudhary et al., 2019); the details of evaluation are provided in the supplementary material. We shall later show (in §3.4.4) that the higher reflections of the velocity field (i.e. ${}_{(2)}\boldsymbol{v}_{(0)}$, ${}_{(3)}\boldsymbol{v}_{(0)},\cdots$) are dispensable for the current problem.

The unknown translational (${{U}_{s\,x}}_{(0)}$) and rotational (${{\Omega}_{s\,y}}_{(0)}$) velocities can be evaluated using the force-free and torque free conditions (24) for first three reflections333Since the second reflection does not contribute to the hydrodynamic drag and torque (Happel & Brenner, 2012; Kim & Karrila, 2013), we evaluate the third reflection to estimate the wall corrections.:

[TABLE]

The $O(\kappa^{2})$ and $O(\kappa^{3})$ corrections are contributions from the wall reflections of the disturbance velocity, which arises when the third reflection is taken into account. The details of the derivation of (37) can be found in the supplementary material. Substitution of (37) in (36) yields $A_{1}\sim O(\kappa^{2}),\;B_{1}\sim H\!a\zeta_{p}O(1),\;\mbox{and\ }C_{1}\sim O(\kappa^{3})$. Since $\kappa\ll 1$, we neglect the $O(\kappa^{2})$ and $O(\kappa^{3})$ wall corrections, and thus use

[TABLE]

This yields the background flow velocity in the frame of reference of $\boldsymbol{U}_{s\,(0)}$ (13)

[TABLE]

The above estimates of $A_{1}$, $C_{1}$ and $U_{s\,x(0)}$ in (38) when substituted in (34) provide an insight into leading order velocity disturbances: (i) stresslet ($\sim O(1/r^{2})$ arising from the resistance to strain), (ii) source-dipole field ($\sim O(1/r^{3})$ arising from the electrophoretic slip), and (iii) octupole field ($\sim O(1/r^{4})$ also arising from the resistance to strain). The other disturbances (multiplying $F_{1},G_{1},H_{1}$), arising due to curvature in the background flow are of $O(\kappa)$ and are accounted in the evaluation of (30).

3.4.4 Evaluation of volume integral

To estimate the migration velocity using (30), we substitute the velocity field $\boldsymbol{v}_{(0)}$ into the integrand. We analyze the intregral for computational convenience and divide the fluid domain ($V_{f}$) into two asymptotic sub-domains (Ho & Leal, 1974, 1976):

[TABLE]

Here, the intermediate radius $\rho$ satisfies: $1\ll\rho\ll 1/\kappa$, $\tilde{r}$ represents the outer coordinates ($\tilde{r}=r\kappa$). The characteristic length scale corresponding to the inner domain ($V_{1}$) is particle radius ‘$a^{*}$’. Whereas, channel width ‘$l^{*}$’ is the characteristic length corresponding to the outer domain ($V_{2}$). Below, we inspect the nature of integrand in each sub-domain and show that the integration in $V_{2}$ is not required to attain the leading order migration velocity.

\justify

Integration in sub-domain $V_{1}$: Substituting (34) and (39) in (19) and performing the integration (30) in $V_{1}$ (which is essentially a spherical shell, provided the particle is not too close to the wall i.e. $s\gg\kappa$, where $\kappa\ll 1$)444The analytical integration of the volume integral is performed using Mathematica 12. The code is provided in the supplementary material., we obtain

[TABLE]

The details of $o(\kappa)$ error incurred while approximating the upper limit of integral to $r=\infty$ is detailed in Appendix C.

In view of the definition of $\beta$ and $\gamma$, the first term in (41) is $O(\kappa)$, and is identical to that obtained by Ho & Leal (1976), who reported the viscoelastic migration of a particle in the absence of electrokinetic effects. The second term is $O(1)$ and is the dominant contribution to the migration in the current study. It is interesting to note that the first term arises purely from the interaction of disturbance field and the bulk flow (i.e. $\boldsymbol{\mathsfbi{w}}^{(1)}_{(0)}\mbox{ and\ }\boldsymbol{\mathsfbi{w}}^{(2)}_{(0)}$). However, the second term arises from the interactions among various disturbances (i.e. $\boldsymbol{\mathsfbi{e}}^{(1)}_{(0)}\bcdot\boldsymbol{\mathsfbi{e}}^{(1)}_{(0)}\mbox{ and\ }\boldsymbol{\mathsfbi{e}}^{(2)}_{(0)}$). This insight and the proportionality to $\beta H\!a\zeta_{p}$ (second term) suggests that the $O(1)$ modification to migration arises from the interaction of stresslet and octupole disturbances (both proportional to $\beta$) with the electrophoresis induced source dipole disturbance (proportional to $H\!a\zeta_{p}$).

\justify

*Integration in sub-domain $V_{2}$: * We now find the contribution of the integral in the outer subdomain $V_{2}$. Since the characteristic length in the outer domain is $l^{*}$, the disturbances are stretched into the outer coordinates ($\tilde{r}=r\kappa$). In the outer coordinate representation, the order of magnitude of the leading order disturbance and test field velocity is

[TABLE]

In the outer coordinates, the undisturbed background flow velocity is

[TABLE]

We substitute (42a) and (43) in (19) and obtain $\tilde{\boldsymbol{\mathsfbi{s}}}_{(0)}$.

The volume integral (30), in $V_{2}$ subdomain, is represented in outer coordinates as

[TABLE]

Upon simplification, we find that the integrand is $O(\kappa^{2})$. Since $\kappa\ll 1$, the dominant contribution can be concluded to arise from the inner integral (i.e. $O(1)$) and the contribution from the outer integral can be neglected. In Appendix D.1, the expression for leading order migration velocity (eq. 41) is verified using an independent calculation following the theorem of Giesekus (1963) for the special case of polymeric fluids ($\Psi_{2}^{*}=-0.5\,\Psi_{1}^{*}$).

4 Results and discussion

The expression for the migration velocity (41) depicts that the direction of lateral migration is determined by the magnitude and sign of the viscometric parameter $\delta\mbox{ and\ }De$. For most viscoelastic fluids $De$ is positive, and $\delta$ is reported to be between -0.5 and -0.7 (Caswell & Schwarz, 1962; Leal, 1975; Koch & Subramanian, 2006; Bird et al., 1987). Since $\gamma$ is negative, the first term in (41) is strictly positive below the centerline ($\beta>0$) and negative above it ($\beta<0$) i.e. the pure elastic stresses push the particle towards the channel axis at $O(\kappa)$. The second term is $O(1)$; its proportionality to $H\!a\zeta_{p}$ indicates that this dominant contribution originates due to the electrophoretic slip i.e., incorporation of electrophoresis enhances the migration by $O(\kappa^{-1})$. For $-0.7<\delta<-0.5$, a particle with positive $\zeta_{p}$ will migrate towards the walls (at $O(1)$), when the electric field is applied in the direction of the flow. The reverse would hold for a particle with negative $\zeta_{p}$ or if the electric field is applied in the opposite direction. These findings qualitatively agree with the experimental observations depicted in Fig.1.

The enhancement demonstrates that an addition of electrokinetic slip on the particle surface can drastically alter the normal stresses which consequently enhances the particle migration. The slip creates a source-dipole disturbance in the fluid which decays rapidly away from the particle ($\sim O(1/r^{3})$). This suggests that the migration occurs from short-ranged interactions, which are unaffected by the boundaries unless the particle is very close to the walls. Furthermore, (41) reveals that the parameters such as wall zeta potential ($\zeta_{w}$) and curvature of the background flow ($\gamma$) do not affect the migration at the present order of approximation.

Since the experimental parameters are in the regime of the theoretical analysis ($De<1$, $De\gg Re_{p}$, $H\!a\zeta_{p}\sim O(1)$, $\kappa\ll 1$), we now quantitatively compare the two. Theoretically, the trajectories of a single particle are calculated using (41) and are converted to an equivalent particle stream-width using $|l^{*}-2(l^{*}-d^{*})|$, for particle focusing at the centerline ($l^{*}$ being the channel height and $d^{*}$ is the particle to wall distance)555In Appendix D.2, the theoretical calculation of trajectory is detailed and are validated with that reported by Ho & Leal (1976) for $H\!a=0$.. We use the relaxation time scale ($(\Psi_{1}^{*}+\Psi_{2}^{*})/\mu^{*}$) corresponding to $C^{*}_{PEO}=250$ ppm as $2.75$ ms (Lu & Xuan, 2015). Figure 3 shows the comparison of particle stream widths for two cases: (i) No electric field (i.e., viscoelastic focusing) and (ii) parallel electric field (i.e., electro-viscoelastic focusing). For case-(i), the theoretical prediction shows modest particle migration. The experiments show negligible focusing as some particles get trapped near the channel corners, which is one of the limitations of viscoelastic focusing at low flowrates (Lu et al., 2017; Tian et al., 2019). Since our model is 2D, it does not account for the trapping of particles at 3D corners. For case-(ii), the theoretical predictions and experimental results confirm an order of magnitude enhancement in focusing. The trapping of particles at the corners is not seen here as the repulsive Maxwell lift is strongest near the corners (Liang et al., 2010). Thus, the agreement between the experiments and theoretical analysis is improved.

Our mathematical framework reveals that the migration in Fig.1(b,c) arises from: (i) Maxwell stresses (eq.26), (ii) pure viscoelastic stresses (eq.41-term1), and (iii) electrophoresis enhanced viscoelastic stresses (eq.41-term2). A comparison of these migration velocities, namely, Maxwell migration (MM), viscoelastic migration (VM) and electro-viscoelastic migration (EVM) is shown in Fig.4. The curves correspond to the case in Fig.1(b): a particle with negative zeta potential and electric field applied parallel to the flow. The migration due to decoupled Maxwell stresses (MM) is dominant near the walls and decays rapidly away from it. Although this wall-repulsive force is negligible in the bulk of the channel (for $\kappa\ll 1$), it helps eliminate the corner equilibrium positions (Liang et al., 2010). The other components VM and EVM act towards the center as depicted by the negative gradient at the central equilibrium position. For conditions corresponding to Fig.1(c), where flow and electric field are opposite to each other ($H\!a<0$), the EVM curve would reverse, as we can see from (41). As a consequence, the central equilibrium position would become unstable, resulting in particle migration towards the walls.

5 Conclusions

Our experimental and theoretical analysis depict an enhancement of $O(\kappa^{-1})$ in migration of a particle in viscoelastic flow, when subjected to a parallel electric field. Here $\kappa$ is the ratio of particle radius to channel size. Using a second-order fluid model, the theoretical analysis reveals the contribution from different underlying stresses around the particle (see eq.26 and eq.41). These constitute an approximate expression for the migration velocity which agrees well with the experimental results. The analytical results further reveal that the particle zeta potential and background shear determine the equilibrium positions. Wall zeta potential, Maxwell stresses, and curvature of the background flow do not affect the migration at the leading order.

Karnis & Mason (1966); Gauthier et al. (1971) suggested that the normal stresses in the undisturbed flow (i.e. tensioned streamlines) are responsible for the migration of the neutral sphere towards the channel center. Later, Ho & Leal (1976) reported that this explanation might not be entirely accurate because of the presence of strong disturbance flow around the sphere. However, the disturbance due to a freely suspended particle in the shear flow does not exhibit a stagnation point, and thus the profile is still predominantly a shear flow (Ho & Leal, 1976; Guazzelli & Morris, 2011). Thus the notion of streamlined tension is not altogether lost and can be applied for a neutral particle. For the case of an electrophoretic particle, the above explanation may not hold as the characteristic disturbance field (being a source dipole) exhibits stagnation points, which should be responsible for the direction of migration addressed in this work. This intriguing feature entails a detailed, perhaps numerical, description of pressure and stress fields around the electrophoretic particle to test this speculation, and warrants further investigation.

The findings of the current work have a direct practical relevance to applications related to Lab-on-a-chip devices used to separate cells in biological fluids which are known to exhibit non-Newtonian behavior (Yuan et al., 2018; Stoecklein & Di Carlo, 2018). The current work addresses the particle migration in dilute polymer concentration regime. For higher polymer concentrations, the direction of migration is reported to be opposite to that of the current work and the explanation for this is currently unknown (Li & Xuan, 2018). Several studies (Rodd et al., 2005; Yang et al., 2011; Kim et al., 2012) have reported that flows consisting semi-dilute polymer concentration regime may exhibit shear-thinning behavior. In such cases, constitutive models like FENE-P and Giesekus fluid can be implemented numerically to accurately predict the particle migration. In our future studies we intend to explore these areas.

The financial support from Indian Ministry of Human Resource Development is gratefully acknowledged.

**Declaration of Interests

**The authors report no conflict of interest.

Appendix A

A.1 Evaluation of $\psi$ using the method of reflections

To satisfy the boundary conditions at both the particle and wall surface, we seek the disturbance potential as a sum of reflections:

[TABLE]

Here, ${}_{(i)}\psi$ represents the ith reflection. The odd reflections satisfy boundary conditions at the particle surface and even reflections satisfy wall boundary conditions. Substituting the above equation in Laplace equation ($\nabla^{2}{\psi}=0$), we obtain for the first reflection:

[TABLE]

and for the second reflection:

[TABLE]

Solution to ${{}_{(1)}}\psi$: Equation (46) suggests that ${{}_{(1)}}\psi$ is a harmonic function which decays as $r\to\infty$. The boundary condition suggests that it must be linear in the ‘driving force’: $\nabla{\phi_{\infty}}(=-{\boldsymbol{e}_{x}})$. Therefore, the solution is sought in terms of spherical solid harmonics (Guazzelli & Morris, 2011) as:

[TABLE]

Solution to ${{}_{(2)}}\psi$: To evaluate ${{}_{(2)}}\psi$, we adopt the approach devised by Faxén (1922) in the context of bounded viscous flows. Using this, the disturbance around the particle is expressed in an integral form which satisfies the boundary conditions at the walls. The first reflection is characterized by particle scale ($a^{*}$). The second reflection ${{}_{(2)}}\psi$ is characterized by a length scale $l^{*}$. Therefore, the coordinates for second reflection are stretched, and are termed as ‘outer’ coordinates. These outer coordinates (denoted by tilde) are defined as:

[TABLE]

Since ${{}_{(2)}}\psi$ depends on ${{}_{(1)}}\psi$ through the wall boundary condition in (47), both reflections must be represented in the coordinates of same scale. The first reflection of potential disturbance in the outer coordinates (${{}_{(1)}}\tilde{\psi}$) is represented as:

[TABLE]

Faxén (1922) represented the fundamental solution of Laplace equation (in Cartesian space) in the terms of Fourier integrals as:

[TABLE]

Here, $\Theta=(\xi\tilde{x}+\eta\tilde{y})/2$ and $\lambda=(\xi^{2}+\eta^{2})^{1/2}$ where $\xi$ and $\eta$ are the variables in Fourier space. We transform the disturbance field (50) by taking derivatives of the above equation (51). For example: (50) contains $\tilde{x}/\tilde{r}^{3}$, which is expressed using the above transformation as:

[TABLE]

This yields an integral representation for the first reflection:

[TABLE]

Here, $b_{1}=1/8$. The boundary condition expressed in (47) suggests a form of solution for ${{}_{(2)}}\tilde{\psi}$ similar to the above equation:

[TABLE]

Substituting (53) and (54) in the boundary condition (47), yields $b_{2}$ and $b_{3}$ in terms of $b_{1}$:

[TABLE]

The next reflection (${{}_{(3)}}\psi$) is $O(\kappa^{3})$. As we restrict $\kappa\ll 1$, the first two reflections capture the leading order contribution to the disturbance potential.

A.2 Evaluation of Maxwell migration

The electric force and torque are expressed as a function of reflections as:

[TABLE]

[TABLE]

The gradient of first reflection ($\bnabla{{}_{(1)}\psi}$) at the particle surface is:

[TABLE]

Here, $\varphi$ and $\theta$ denote the azimuthal and polar angle, respectively. As the wall reflection of the potential is represented in the outer scaled coordinates, $\bnabla{{}_{(2)}\psi}$ at the particle surface ($r=1$) is expressed as a Taylor series expansion of $\kappa\tilde{\bnabla}\tilde{\psi_{2}}$ about the origin ($\tilde{r}\rightarrow 0$). Here $\tilde{\bnabla}$ is the spatial gradient in outer scaled coordinates: $\left\{{\frac{\partial}{\partial\tilde{x}},\,\frac{\partial}{\partial\tilde{y}},\,\frac{\partial}{\partial\tilde{z}}}\right\}$. Since the particle is small with respect to the wall-particle gap ($\kappa/s\ll 1$), we retain only the zeroth and first order terms:

[TABLE]

Substituting ${{}_{(1)}\psi}$ $\&$ ${{}_{(2)}\psi}$ in both (56) $\&$ (57) and using (58) $\&$ (59), we obtain:

[TABLE]

Appendix B Test field ($\boldsymbol{u}^{t}$)

The test field is chosen to be that generated by a sphere moving in the positive z-direction with a unit velocity in a quiescent Newtonian medium.

[TABLE]

To solve for $\boldsymbol{u}^{t}$, we again make use of the fact that the integration in $V_{1}$ is sufficient for the leading order solution and solve for the first reflection of $\boldsymbol{u}^{t}$ using Lamb’s solution:

[TABLE]

Appendix C

Here we analyze the integrand of (30) in sub-domain ${V_{1}}=\left\{{\boldsymbol{r}:1\leq|r|\leq\rho}\right\}$, where $\rho$ is an intermediate radius ($1\ll\rho\ll\kappa^{-1}$). The order of magnitude of the leading order disturbance and test field velocity are:

[TABLE]

From (39), we obtain the order of magnitude of the undisturbed velocity

[TABLE]

We substitute (63 & 64) in (19) and use it in (30) to find the order of magnitude of integrand $\boldsymbol{\mathsfbi{s}}_{(0)}\colon\bnabla\boldsymbol{u}^{t}$ in $V_{1}$:

[TABLE]

The coefficients of different terms in (65) depict contributions to the integrand arising from various interactions. The radial order describes the range of interaction. Equation (41) suggests that the contribution from all the terms in (65) vanish except the second and fifth term.

We now analyze the error associated with the approximation of the upper limit of integral as $r=\infty$ in (41). Following Ho & Leal (1974), the intermediate radius $\rho$ is written as

[TABLE]

Out of the two non-vanishing components of integrand in (65), the fifth term decays the slowest $\sim O(1/r^{4})$. Thus, the leading order error is generated due to it. Its contribution to the integral ($\int_{V_{1}}\boldsymbol{\mathsfbi{s}}_{0}:\bnabla{\boldsymbol{u}^{t}}dV$) is

[TABLE]

Using the definitions in (15), we conclude that the error incurred above is $O(\kappa^{2-\eta})$. For $0<\eta<1$ and $1\ll\rho\ll\kappa^{-1}$, the error is $o(\kappa)$.

Appendix D

D.1 Validation using the 3-D Giesekus theorem

For a special case of $\delta=-1$ (corresponding to $\Psi_{2}^{*}=-0.5\,\Psi_{1}^{*}$), Giesekus theorem (Giesekus, 1963; Bird et al., 1987) states that the velocity in a second-order fluid is unaltered from the Newtonian velocity. However, the steady pressure field alters as

[TABLE]

We substitute the Newtonian velocity and pressure field and evaluate the modified pressure. This pressure is now substituted in the total hydrodynamic stress tensor $\boldsymbol{\mathsfbi{T}}=-P\boldsymbol{\mathsfbi{I}}+\boldsymbol{\mathsfbi{E}}^{(1)}+De\,\boldsymbol{\mathsfbi{S}}.$ The polymeric stress ($\boldsymbol{\mathsfbi{S}}$) can be represented as an addition of the co-rotational derivative ($\boldsymbol{\mathsfbi{S}}_{C}$) and quadratic component ($\boldsymbol{\mathsfbi{S}}_{Q}$) of the rate of strain tensor (Koch & Subramanian, 2006)

[TABLE]

where $\boldsymbol{\mathsfbi{R}}^{(1)}$ is the rate of rotation tensor $(\bnabla\boldsymbol{U}-\bnabla\boldsymbol{U}^{\dagger})$. In the limit of $\delta=-1$, the quadratic stress vanishes.

Accounting for the leading order effects of viscoelasticity on the lateral force (provided $\delta=-1$) yields

[TABLE]

Here, $\boldsymbol{n}$ is the unit normal vector pointing outside the sphere and steady-state $\boldsymbol{\mathsfbi{S}}_{C\,(0)}$ is evaluated by substituting the Newtonian field ($\boldsymbol{U}_{(0)},P_{(0)}$) in (69). Solving the above integral results in the following lateral force

[TABLE]

which is identical to that obtained from (41) for $\delta=-1$. In the above calculation, $De\,\delta(-\frac{24\pi}{3}\beta\gamma-12\pi\beta H\!a\zeta_{p})$ is the contribution from pressure field, and $De\,\delta(\frac{44\pi}{3}\beta\gamma+12\pi\beta H\!a\zeta_{p})$ is the contribution from the co-rotational stress. The analysis here shows that the enhancement in migration comes entirely from the quadratic stress.

D.2 Evaluation and validation of the particle trajectory

The ratio of migration velocity to the translational velocity is

[TABLE]

At the leading order, $U_{s\,x}$ can be approximated as $U_{s\,x(0)}$ because the corrections arrive at $O(De)$ (due to weak non-Newtonian effects) and $O(\kappa^{2})$ (due to wall induced viscous resistance). $U_{mig}$ has both hydrodynamic and Maxwell stress contributions (i.e. $U_{mig}^{H}+U_{mig}^{M}$). Using first order Euler method and substituting $U_{mig}^{H}$, $U_{mig}^{M}$, and $U_{s\,x(0)}$, we obtain the equation governing the particle trajectory

[TABLE]

Since $\alpha\sim O(\kappa^{-1})$, $\gamma\sim O(\kappa)$, $H\!a(\zeta_{p}-\zeta_{w})\sim O(1)$, the trajectory equation can be approximately represented as

[TABLE]

We convert this trajectory into an equivalent particle stream width (PSW) for the case of centerline migration (fig.1-b). Since the concentration of suspension is dilute ($<$0.1$\%$ in volume fraction), we neglect particle-particle interactions. At the inlet, particles are homogeneously suspended and cover the entire width of the channel. As the suspension flows, the particles migrate towards the centerline.

Experimentally, PSW is measured by tracking the spread of particles (i.e. for centerline focusing, the particle closest to the walls is tracked). Theoretically, using (74), we find the trajectory of this particle, which is initially suspended near the wall ($s^{1}=0.95$). The non-dimensional PSW associated with this trajectory is evaluated using $|1-2(1-s^{n})|$. For instance, $s^{n}=0.8\mbox{ or\ }0.2$ corresponds to PSW$=0.6$ i.e. the dilute suspension of particles cover 60% of the channel width. In dimensional form, the particle stream width is $|l^{*}-2(l^{*}-d^{*\,n})|$. Fig. 5 depicts the stream width of particles suspended in a pressure driven flow of second-order fluid, in the absence of an electric field. The result matches well with that reported by (Ho & Leal, 1976, p.796) (the trajectories are converted to PSW using $|l^{*}-2(l^{*}-d^{*\,n})|$). The dashed line shows the particle stream width when wall effects are taken into account (at the leading order). Results in §4 include this first order wall correction.