The motion of a neutrally buoyant particle of an elliptic shape in two dimensional shear flow: a numerical study

Shih-Lin Huang, Shih-Di Chen, Tsorng-Whay Pan^{1}^{1}1Corresponding authors. e-mail: pan@math.uh.edu, mechang@mail.iam.ntu.edu.tw,
Chien-Cheng Chang, Chin-Chou Chu Institute of Applied Mechanics, National Taiwan University, Taipei 106, Taiwan, Republic of China Center for Advanced Studies in Theoretical Sciences, National Taiwan University, Taipei 106, Taiwan, Republic of ChinaDepartment of Mathematics, University of Houston, Houston, Texas 77204, USA

Abstract In this paper, we investigate the motion of a neutrally buoyant cylinder of an elliptic shape freely moving in two dimensional shear flow by direct numerical simulation. An elliptic shape cylinder in shear flow, when initially being placed at the middle between two walls, either keeps rotating or has a stationary inclination angle depending on the particle Reynolds number , where is the shear rate, is the semi-long axis of the elliptic cylinder and is the kinetic viscosity of the fluid. The critical particle Reynolds number for the transition from a rotating motion to a stationary orientation depends on the aspect ratio and the confined ratio where is the semi-short axis of the elliptic cylinder and is the distance between two walls. Although the increasing of either parameters makes an increase in , the dynamic mechanism is distinct. The variation causes the change of geometry shape; however, the variation influences the wall effect. The stationary inclination angle of non-rotating slender elliptic cylinder with smaller confined ratio seems to depend only on the value of . An expected equilibrium position of the cylinder mass center in shear flow is the centerline between two walls, but when placing the particle away from the centerline initially, it migrates either toward an equilibrium height away from the middle between two walls or back to the middle depending on the confined ratio and particle Reynolds number. keywords: Shear flow; neutrally buoyant elliptic cylinder; equilibrium height; stationary inclination angle; fictitious domain/distributed Lagrange multiplier method.

## 1 Introduction

The problem of particle motion in shear flow is crucially important in many engineering fields such as the handling of a fluid-solid mixture in slurry, colloid, and fluidized bed. Experimentally Segré and Silberberg [22, 23] studied the migration of neutrally buoyant spheres in a tube Poiseuille flow and obtained that the particles migrate away from the wall and the centerline to accumulate at about 0.6 of the tube radius from the central axis. The experiments of Segré and Silberberg [22, 23] have had a large influence on fluid mechanics studies of migration and lift of particles. Comprehensive reviews of experimental and theoretical works have been given by Brenner [1], Cox and Mason [6], Feuillebois [10] and Leal [15].

Jeffery [14] analyzed the equation of motion of a particle immersed in an unbounded viscous fluid under Stokes flow, where the fluid and particle inertia could be completely negligible relative to viscous forces. According to the equation of motion, he corroborated a periodic rotation of an ellipsoid in a simple shear flow. Concerning the theoretical results of the neutrally buoyant particle migration in linear shear flow, Bretherton [2] found an expression for the lift force per unit length on a cylinder in an unbounded two-dimensional linear shear flow at small Reynolds number. Saffman’s lift force [21] on a sphere of radius in an unbounded linear shear flow with shear rate is where is the kinetic viscosity of the fluid, is the density of the fluid, is the particle Reynolds number, and is the slip velocity of the sphere. In a bounded linear shear flow, Ho and Leal [13] examined the motion of a rigid sphere with inclusion of the inertia effects at small Reynolds numbers by a regular perturbation method. The sphere reaches a stable lateral equilibrium height which is the midway between the walls. Vasseur and Cox [24] also obtained the same stable lateral equilibrium height. Ho and Leal require that which is more restrictive than the one required by Vasseur and Cox where is the confined ratio, being the distance between two walls. Via direct numerical simulation, Feng et al. [9] investigated the motion of neutrally buoyant and non-neutrally buoyant circular particle in plane shear and Poiseuille flows using a finite element method and obtained qualitative agreement with the results of perturbation theories and of experiments. The numerical results of a neutrally buoyant circular cylinder in a shear flow of have been discussed in details. The cylinder migrates back to the midway between two walls due to the wall repulsion at the small Reynolds number. They have suggested that three factors, namely the wall repulsion due to a lubrication effect, the slip velocity, and the Magnus type of lift, are possibly responsible for the lateral migration. Ding and Aidun [7] studied numerically the dynamics of a cylinder of circular or elliptic shape suspended in shear flow at various particle Reynolds number. They obtained the transient from being rotary to stationary as the particle Reynolds number is increased for an elliptic cylinder. Zettner and Yoda [26] followed up the work done by Ding and Aidun [7] and tested the neutrally buoyant cylinders of elliptic and non-elliptic shape over a wide range of aspect ratios within a moderate range of Re and obtained the rotational motion and the stationary orientation behavior of a neutrally buoyant rigid body in a simple shear flow qualitatively agree with the results of numerical simulations in a substance.

In this paper, we have applied a distributed Lagrange multiplier/fictitious domain method developed in [4, 20] to study the dynamics of a neutrally buoyant elliptic cylinder in two-dimensional shear flow. Such methodology has been validated numerically in [4, 20] by comparing with the computational results by Ding and Aidun in [7] for a cylinder of circular and elliptic shape and the experimental results by Zettner and Yoda in [27] for a circular cylinder. The dynamics of a neutrally buoyant particle of elliptic shape has been studied here based on three physical parameters: the aspect ratio of the particle , the degree of confinement , and the particle Reynolds number . Then we investigate on the equilibrium height of a neutrally buoyant elliptic cylinder in shear flow. Feng and Michaelides [8] have investigated the equilibrium heights of non-neutrally and almost neutrally buoyant circular cylinders in two-dimensional shear flow. In their simulations, the density ratio between the solid and fluid is between 1.005 and 1.1. The equilibrium heights of their lightest circular cylinder (the density ratio of 1.005) are far below the centerline for between 2 and 4.5. In [20], Pan et al. have obtained equilibrium positions off the centerline between two walls, besides the expected one in the middle between two walls, for a neutrally buoyant circular cylinder at different Reynolds numbers in shear flow. Such new equilibrium height depends on the particle Reynolds number and the confined ratio . For a neutrally buoyant cylinder of elliptic shape in shear flow, we have obtained not just similar results when the initial position of the mass center is off the middle between two walls. At higher particle Reynolds number, it is not surprised to obtain that the neutrally buoyant cylinder of elliptic shape at the off-the-middle equilibrium height can keep a stationary inclination angle as those in the middle between two walls since it is still in a shear flow.

The content of this paper is as follows: In Section 2 we briefly introduce a fictitious domain formulations of the model problem and computational methods associated with the neutrally buoyant long particle cases; then in Section 3 we present and discuss the numerical results. The conclusions are summarized in Section 4.

## 2 A fictitious domain formulation of the model problem

A fictitious domain formulation with distributed Lagrange multipliers for flow around freely moving particles and its associated computational methods have been developed and tested in, e.g., [11, 12, 16, 17, 18, 19, 25]. For the cases of a neutrally buoyant elliptic cylinder in two dimensional flows, the methodology has been developed and validated in [4, 20]. Let be a rectangular region filled with a Newtonian viscous incompressible fluid (of density and dynamic viscosity ) and containing a freely moving neutrally buoyant rigid particle centered at of density , as depicted in Figure 1. The flow is modeled by the Navier–Stokes equations and the motion of the particle is modeled by the Euler–Newton’s equations. The basic idea of the fictitious domain method is to imagine that the fluid fills the entire space inside as well as outside the particle boundary. The fluid–flow problem is then posed on a larger domain (the “fictitious domain”). The fluid inside the particle boundary must exhibit a rigid–body motion. This constraint is enforced using the distributed Lagrange multiplier, which represents the additional body force per unit volume needed to maintain the rigid–body motion inside the particle boundary, much like the pressure in incompressible fluid flow, whose gradient is the force required to maintain the constraint of incompressibility. For flow around a freely moving neutrally buoyant particle of a long body shape, the fictitious domain formulation with distributed Lagrange multipliers is as follows

(1) | |||

(2) | |||

(3) | |||

(4) | |||

(5) | |||

(6) | |||

(7) |

where and denote velocity and pressure, respectively, is a Lagrange multiplier, is the translation velocity of the particle , is the angular velocity of , and is the angle between the horizontal direction and the long axis of the elliptic cylinder (see Fig. 1). We suppose that the no-slip condition holds on . We also use, if necessary, the notation for the function . The function spaces in equations (1)-(7) are defined by

with , , and an inner product on which can be the standard inner product on . For simple shear flow, we have on the top wall and on the bottom wall.

###### Remark 1.

The hydrodynamical forces and torque imposed on the rigid body by the fluid are built in equations (1)-(7) implicitly (see [11, 12] for details), thus we do not need to compute them explicitly in the simulation. Since in equations (1)-(7) the flow field is defined on the entire domain , it can be computed with a simple structured grid.

###### Remark 2.

###### Remark 3.

The method of numerical solution is actually a combination of a distributed Lagrange multiplier based fictitious domain method and an operator splitting method. For space discretization, we use –iso– and finite elements for the velocity field and pressure, respectively (like in [3]). In time advancing, we apply the Lie’s scheme in [5] to obtain a sequence of sub–problems for each time step (see [17, 20] for details). The computational method has been validated in [4, 20] by comparing with the computational results by Ding and Aidun in [7] for a cylinder of circular and elliptic shape and the experimental results by Zettner and Yoda in [27] for a circular cylinder.

## 3 Results and discussions

### 3.1 A neutrally buoyant elliptic cylinder placed initially at the middle between two walls

The motion of a neutrally buoyant elliptic cylinder placed initially at the middle between two walls in two–dimensional shear flow was studied by Ding and Aidun in [7] via a lattice–Boltzmann method. When the particle Reynolds number is increasing to and less than the critical value, , the elliptic shape particle keeps rotating at the middle between two walls but the period of the rotation becomes longer and longer. Once is beyond the critical value, the elliptic shape cylinder stops rotating but has a stationary inclination angle with respect to the wall direction since the torques before and after the cylinder are strong enough to hold the particle of long body shape (e.g., see Fig. 2 and [7]). Such orientation in shear flow with a stationary inclination angle is a quite different behavior. We have studied the dynamics of such motion based on the particle Reynolds number , the confined ration , and the aspect ratio in the following.

In a computational domain , a neutrally buoyant elliptic cylinder is placed at the middle between two walls initially. Thus the distance between two walls is . The fluid density is . The shear rate is fixed at so for a given particle Reynolds number, the kinetic viscosity is given by . At first, we focus on the rotation of an elliptic cylinder as . Fig. 3 shows the relation between the dimensionless period and the Reynolds number . The effect on the long body with a fixed in Fig. 3 shows that a slender ellipse reveals a more rapid increase of within a small Re than a fat ellipse does within a large , and its is much smaller than the one associated with a fat ellipse, which is closer to a circular shape. In addition, in Fig. 3, regarding the effect under a fixed aspect ratio and 0.5, results show that the less confined long body also present a more rapid increase of than those of the more confined ones and exhibit a smaller . The minimum angular velocity and the Reynolds number do show a linear relation as in Fig. 4 for the Reynolds number slightly less than the critical value. Using such relation, the critical Reynolds number is predicted here and in [7]. Although the increasing of either parameters makes an increase in (e.g., see Figs. 3 and 4), the dynamic mechanism is distinct. The variation causes the change of geometry shape; however, the variation influences the wall effect.

Concerning the angular velocity of an elliptic cylinder suspended in shear flow, the minimal angular velocity decreases rapidly to about zero when the particle Reynolds number is closer but less than , e.g., as in Fig. 5, since the motion of the elliptic cylinder is about to transit into the one with a stationary orientation in shear flow. But its maximal angular velocity shown in Fig. 5 have different behavior since the maximal angular velocity happens when the direction of the long axis is about perpendicular to the shear direction. For the small particle Reynolds numbers, the maximal and minimal values of the angular velocity are very close to the Jeffery’s solution as in Fig. 5.

As , an elliptic cylinder in shear flow has a stationary inclination angle. We have obtained that the inclination angle seems to depend only on the value of for the aspect ratio less than or equal to 0.5 and the confined ratio less than or equal to 0.4 as shown in Fig. 6. But for a fat elliptic cylinder (i.e., ) or an elliptic cylinder with higher confined ratio , the inclination angle is smaller. The increase in stationary inclination angle with increasing have also been obtained by Ding and Aidun in [7] and Zettner and Yoda in [27]. Qualitatively, we have obtained similar behavior.

### 3.2 A neutrally buoyant elliptic cylinder placed initially away from the midway between two walls

In Ho and Leal [13] and Vasseur and Cox [24], they concluded that the sphere reaches a stable lateral equilibrium position which is the midway between the walls for small particle Reynolds numbers. In Feng et al. [9], a circular cylinder migrates back to the midway between two walls at when placing it away from the middle between two walls. Feng et al. have suggested that three factors, namely the wall repulsion due to a wall repulsion force, the slip velocity, and the Magnus type of lift, are possibly responsible for the lateral migration. Feng and Michaelides [8] have investigated the equilibrium heights of non-neutrally buoyant circular cylinders in two-dimensional shear flow. In their simulations, the density ratio between the solid and fluid is between 1.005 and 1.1. The equilibrium heights of their lighter circular cylinder (the density ratio of 1.005) are far below the centerline for between 2 and 4.5. In [20], Pan et al. have obtained another equilibrium height which is off the middle between two walls for a neutrall buoyant circular cylinder in shear flow, besides the one in the middle between two walls. Such off-the-middle equilibrium height does depend on the particle Reynolds number and the confined ratio . In the previous section, we have obtained that the midway is always (at least in the range we have studied in this paper) the equilibrium height for the cylinder of elliptic shape when it is placed there initially. For the study of the off-the-middle equilibrium height, we have considered a neutrally buoyant elliptic cylinder in with , , and or 1. When placing the mass center initially 0.1 length unit below the middle between two walls, the mass center of the freely moving cylinder migrates back to the middle between two walls for as in Fig. 7; but for higher particle Reynolds numbers, the mass center migrates to an equilibrium height away from the middle as shown in Figs. 7, 8 and 10. It is not surprised to find out that, like those staying at the middle, the cylinder either keeps rotating or has a stationary inclination angle when being away from the midway since it is still in a shear flow, e.g., see Fig. 8. For those away from the middle with a stationary orientation, the inclination angles are 17.62, 20.69 and 22.18 degrees for the particle Reynolds numbers , 13, and 15, respectively. Those angles are three to four degrees smaller than the ones associated with the cylinders at the middle between two walls at the same Reynolds numbers. From the snapshots in Fig. 8 for , we can see that the smaller inclination angle is due to the effect of the confinement from the bottom wall, which is consistent with the results discussed in the previous section.

About the angular velocity, both cylinders keep rotating at as shown in Fig. 10 and the period of the one off the middle initially is slightly smaller due to the effect of the confinement, which is also consistent with the results shown in Fig. 3. At , the histories of the angular velocity of the cylinder at two different initial positions are almost identical. They both rotate right away to their stationary inclination angles, respectively, and then stop rotating as indicated in Fig. 10. For the one in the midway between two walls, it remains there due to the symmetry with respect to the mass center. But for the other one not in the middle between two walls initially, Figs. 10 and 10 show that the cylinder moves with a stationary inclination angle in shear flow toward the bottom wall slowly and approaches its equilibrium height at which, we believe, the force from the effect of slip velocity and the wall repulsive force are balanced as discussed in the following.

The numerical study on circular cylinders in [20] about the effect of the Magnus lift force was done by considering a circular cylinder without rotation in shear flow, which can be achieved by adding a constraint of zero rotation velocity numerically. Pan et al. obtained that the mass center of a circular cylinder with no rotation is lower than that of the one with totally free motion, which is exactly the effect of the Magnus lift force. In this article, we have studied the Magnus lift force by comparing the height of the elliptic cylinder with a fixed orientation angle to the freely rotating one at . When the inclination angle reaches the specified value after , the ellipse cylinder is not allowed to rotate any more in shear flow just like the one done in [20]. The chosen inclination angles are 10, 15, 20, and 30 degrees. The histories of the height of the mass center in Fig. 11 show that the cylinder without rotation goes to a lower equilibrium height, which indicates that the Magnus lift force from the rotation does lift a neutrally buoyant elliptic shape cylinder in shear flow.

Concerning the effect of slip velocity, the cylinder is moving to the left in the lower region of the computational domain in shear flow (see Fig.1 for the set-up of the boundary conditions) due to the initial position of the cylinder which is below the centerline. For getting the slip velocity, we first compute the fluid horizontal speed on the line through the mass center in front of the cylinder at the distance of the half of the computational domain width and then minus the horizontal speed of the cylinder to obtain the slip velocity. At and 15, the slip velocities of the cylinder are almost zeros for those initially placed at the middle between two walls as in Fig. 10. But for the ones placed initially below the centerline at and 15, both slip velocities become negative after a very short initial transition period. The negative slip velocity means that the cylinder horizontal speed to the left is slower than the fluid speed to the left. Thus the cylinder lags the fluid and then is pushed toward the region with faster flow speed which is the region next to the bottom wall for both and 15.

Overall the results in Fig. 7 show that, for (resp., ), the combined effect of the wall repulsion force and the Magnus lift force is relatively stronger than the effect of the slip velocity so that the cylinder is pushed back to the midway for (resp., ). But for (resp., ), the effect of the slip velocity dominates the other two so that the cylinder stays away from the middle as (resp., ). For those having stationary inclination angle, the balance between the effect of the slip velocity and the wall repulsion force determines the off-the-middle equilibrium height of the cylinder in shear flow when the initial position is not at the middle of two walls. For the actual distance from the cylinder mass center to the bottom wall, Fig. 12 shows that the off-the-middle equilibrium heights are about the same for for both values of . The off-the-middle height is a function of the Reynolds number, for (resp., for ) as (resp., ).

## 4 Conclusions

We have investigated the motion of a neutrally buoyant cylinder of an elliptic shape in two dimensional shear flow of a Newtonian fluid by direct numerical simulation. An elliptic shape cylinder in shear flow, when initially being placed at the middle between two walls, either keeps rotating or has a stationary inclination angle depending on the particle Reynolds number . The critical particle Reynolds number for the transition from a rotating motion to a stationary orientation depends on the aspect ratio and the confined ratio . Although the increasing of either parameters makes an increase in , the dynamic mechanism is distinct. The variation causes the change of geometry shape; however, the variation influences the wall effect. The stationary inclination angle of non-rotating slender elliptic cylinder with smaller confined ratio seems to depend only on the value of . An expected equilibrium position of the cylinder mass center in shear flow is the centerline between two walls; but when placing the particle away from the centerline initially, it either migrates back to the midway or moves away from the middle between two walls depending on the particle Reynolds number and the confined ratio. For those having stationary inclination angle, the balance between the effect of the slip velocity and the wall repulsion force does play a role for determining the equilibrium height of the cylinder in shear flow when the initial position is not at the middle of two walls.

Acknowledgments. T.-W. Pan acknowledges the support by the US NSF under Grant No. DMS-0914788. S.-L. Huang, S.-D. Chen, C.-C. Chu, C.-C. Chang acknowledge the support by the National Science Council (Taiwan, ROC) under Contract Numbers, NSC97-2221-E-002-223-MY3, NSC99-2628-M-002-003 and NSC100-2221-E-002-152-MY3.

## References

- [1] H. Brenner, “Hydrodynamic resistance of particles at small Reynolds numbers,” Adv. Chem. Engng. 6, 287 (1966).
- [2] F. P. Bretherton, “The motion of rigid particles in a shear flow at low Reynolds number,” J. Fluid Mech. 14, 284 (1962).
- [3] M. O. Bristeau, R. Glowinski, J. Periaux, “Numerical methods for the Navier-Stokes equations. Applications to the simulation of compressible and incompressible viscous flow,” Computer Physics Reports 6, 73 (1987).
- [4] S.-D. Chen, T.-W. Pan, C.-C. Chang, “The motion of a single and multiple neutrally buoyant elliptical cylinders in plane Poiseuille flow,” Phys. Fluids 24, 103302 (2012).
- [5] J. Chorin, T. J. R. Hughes, M. F. McCracken, J. E. Marsden, “Product formulas and numerical algorithms,” Commun. Pure Appl. Math. 31, 205 (1978).
- [6] R. G. Cox, S. G. Mason, “Suspended particles in fluid flow through tubes,” Ann. Rev. Fluid Mech. 3, 291 (1971).
- [7] E. Ding, C. K. Aidun, “The dynamics and scaling law for particles suspended in shear flow with inertia,” J. Fluid Mech. 423, 317 (2000).
- [8] Z.-G. Feng, E. E. Michaelides, “Equilibrium position for a particle in a horizontal shear flow,” Int. J. Multiphase Flow 29, 943 (2003).
- [9] J. Feng, H. H. Hu, D. D. Joseph, “Direct simulation of initial value problems for the motion of solid bodies in a Newtonian fluid. Part 2: Couette and Poiseuille flows,” J. Fluid Mech. 277, 271 (1994).
- [10] F. Feuillebois, “Some theoretical results for the motion of solid spherical particles in a viscous fluid”, Multiphase science and technology (edited by G.F. Hewitt, J.M. Delhaye and N. Zuber, Hemisphere Pub. Corp., New York, 1989), Vol. 4, 583-798.
- [11] R. Glowinski, T.-W. Pan, T. Hesla, D. D. Joseph, “A distributed Lagrange multiplier/fictitious domain method for particulate flows,” Int. J. Multiphase Flow 25, 755 (1999).
- [12] R. Glowinski, T.-W. Pan, T. Hesla, D. D. Joseph, J. Periaux, “A fictitious domain approach to the direct numerical simulation of incompressible viscous flow past moving rigid bodies: Application to particulate flow,” J. Comput. Phys. 169, 363 (2001).
- [13] B. P. Ho, L. G. Leal, “Inertial migration of rigid spheres in two-dimensional unidirectional flows,” J. Fluid Mech. 65, 365 (1974).
- [14] G. B. Jeffery, “The motion of ellipsoidal particles immersed in a viscous fluid,” Proc. R. Soc. Lond. A 102, 161 (1922).
- [15] L. G. Leal, “Particle motions in a viscous fluid,” Ann. Rev. Fluid Mech. 12, 435 (1980).
- [16] T.-W. Pan, C.-C. Chang, R. Glowinski, “On the motion of a neutrally buoyant ellipsoid in a three-dimensional Poiseuille flow,” Comput. Methods Appl. Mech. Engrg. 197, 2198 (2008).
- [17] T.-W. Pan, R. Glowinski, “Direct simulation of the motion of neutrally buoyant circular cylinders in plane Poiseuille flow,” J. Comput. Phys. 181, 260 (2002).
- [18] T.-W. Pan, R. Glowinski, “Direct simulation of the motion of neutrally buoyant balls in a three-dimensional Poiseuille flow,” C. R. Mecanique, Acad. Sci. Paris 333, 884 (2005).
- [19] T.-W. Pan, D. D. Joseph, R. Bai, R. Glowinski, V. Sarin, “Fluidization of 1204 spheres: simulation and experiments,” J. Fluid Mech. 451, 169 (2002).
- [20] T.-W. Pan, S.-L. Huang, S.-D. Chen, C.-C. Chu, C.-C. Chang, “A numerical study of the motion of a neitrally buoyant cylinder in two dimensional shear flow,” Computers & Fluids 87, 57 (2013).
- [21] G. Saffman, “The lift on a small sphere in a slow shear flow,” J. Fluid Mech. 22, 385 (1965).
- [22] G. Segré, A. Silberberg, “Radial particle displacements in Poiseuille flow of suspensions,” Nature 189, 209 (1961).
- [23] G. Segré, A. Silberberg, “Behavior of macroscopic rigid spheres in Poiseuille flow. Part I,” J. Fluid Mech. 14, 115 (1962).
- [24] P. Vasseur, R. G. Cox, “The lateral migration of a spherical particle in two-dimensional shear flows,” J . Fluid Mech. 78, 385 (1976).
- [25] B. H. Yang, J. Wang, D. D. Joseph, H. H. Hu, T.-W. Pan, R. Glowinski, “Numerical study of particle migration in tube and plane Poiseuille flows,” J. Fluid Mech. 540, 109 (2005).
- [26] C. M. Zettner, M. Yoda, “Moderate-aspect-ratio elliptical cylinders in simple shear with inertia,” J. Fluid Mech. 442, 241 (2001).
- [27] C. M. Zettner, M. Yoda, “The circular cylinder in simple shear at moderate Reynolds numbers: An experimental study,” Expts. Fluids 30, 346 (2001).