Run Li,Jingsheng Zhang ,2,Yumei Yong ,Yang Wang ,*,Chao Yang ,*
1 Key Laboratory of Green Process and Engineering,Institute of Process Engineering,Chinese Academy of Sciences,Beijing 100190,China
2 Beijing Research Institute of Chemical Industry,SINOPEC,Beijing 100013,China
Keywords:Shear flow Liquid sphere Numerical simulation Streamline Jeffery orbit
ABSTRACT This work presents a numerical investigation on steady internal,external and surface flows of a liquid sphere immersed in a simple shear flow at low and intermediate Reynolds numbers.The control volume formulation is adopted to solve the governing equations of two-phase flow in a 3-D spherical coordinate system.Numerical results show that the streamlines for Re=0 are closed Jeffery orbits on the surface of the liquid sphere,and also closed curves outside and inside the liquid sphere.However,the streamlines have intricate and non-closed structures for Re≠0.The flow structure is dependent on the values of Reynolds number and interior-to-exterior viscosity ratio.
Shear flow generally exists in multiphase dispersions in process industry.In some practical cases,a liquid sphere(droplet)subjects to the force of gravity and shear simultaneously.For example,the droplets in a stirred tank would bear shear force from the continuous phase at low and moderate Re(e.g.,in high viscous systems and fine emulsions).Convection and shear have appreciable effects on the transport processes at moderate Reynolds numbers.A thorough understanding of the flow structure around a single droplet in simple shear flow would help in gaining insight into the transport process of a droplet in pure shear and complex flows.The research on pure shear flow at low and moderate Re would also help us to understand the complex interaction of shear and advection co-existing in general liquid-liquid systems.The internal and external flow fields of dispersed phase particles(including bubble and drop)would display special flow structures in shear flows.Thus it is necessary to study the fluid mechanics of a single particle in shear flows for extensive understanding of rheological properties of multiphase dispersions.On the other hand,mass and heat transfer correlates closely with the inside and outside flow fields.This motivates the present study on the flow structures around single particles.
Peery[1]used a singular perturbation technique to study the effect of weak fluid inertia on the fluid velocity field around a rigid or deformable sphere in simple shear flow.Roberson and Acrivos[2]investigated theoretically and experimentally the fluid velocity field around a freely suspended cylinder in simple shear flow at low Reynolds numbers and found that the region of closed streamlines had a finite extent along the direction of flow.Poe and Acrivos[3]studied a solid sphere rotating freely in simple shear flow experimentally for moderate values of Reynolds number up to 10 and obtained the rotation rates of such a sphere.Subramanian and Koch[4]deduced that the fluid inertia made streamlines near a solid sphere open at non-zero Re,instead of remaining closed as in the case of Re=0.For non-zero but very small Re,centrifugal forces caused the streamlines in the flow-gradient plane spiral away from the particle surface.Mikulencak and Morris[5]quantified the particle rotation rate for a solid sphere in simple shear flow and its contribution to the fluid stress by using a finite element method.Yang et al.[6]studied similar problem numerically by a finite difference method and their numerical results of particle rotation rate were consistent with those of Mikulencak and Morris[5].These theoretical analyses[5,6]were limited to Re?1 and solid spheres.Liquid spheres immersed in simple shear flow at intermediate Re are seldom targeted.The numerical results of velocity and stress[6]have been used to analyze the stresslet for liquid spheres[7].Mao et al.[8]investigated numerically the fully developed steady flow of non-Newtonian yield viscoplastic fluid through concentric and eccentric annuli.The fluid rheology is described with the Herschel-Bulkley model.The numerical simulation based on a continuous viscoplastic approach to the Herschel-Bulkley model is found in poor accordance with the experimental data on volumetric flow rate of a bentonite suspension.Fan and Yin[9]investigate the interaction of two bubbles rising side by side in shear-thinning fluid using volume of fluid(VOF)method coupled with continuous surface force(CSF)method.By considering rheological characteristics of fluid,this approach was able to accurately capture the deformation of bubble interface,and validated by comparing with the experimental results.
In this work,we determine numerically the flow field around a neutrally buoyant liquid sphere in simple shear flow at finite Reynolds numbers with a control volume formulation.We present in detail the flow fields inside and outside a liquid sphere in simple shear and the flow structure on drop surface.It is believed that the flow features revealed by numerical simulation will be useful for further analysis of heat and mass transfer as well as liquid-liquid chemical reactions.
A rigid liquid sphere is placed at the origin of coordinate system in a Newtonian fluid,which is subject to simple constant shear far from the droplet.The continuous and dispersed phases have equal densities.The flow field without the central sphere is given in the Cartesian coordinates as,whereis the velocity gradient of simple shear flow.Flow circulations may exist inside a liquid sphere due to the shear stress from the continuous phase.The physical properties of two phases may be different.In the laminar flow regime,the velocity(u)and pressure(p)in each phase are governed by the continuity and Navier-Stokes equations,in dimensionless form,as follows

where the subscript i=1 is for the droplet and i=2 for the continuous phase.Coordinates are non-dimensionalized by liquid sphere radius a,velocity by a,and stresses bywhereμ is the viscosity of continuous phase.The Reynolds number is defined by,the viscosity ratio by λ = μ1/μ2,and ρ is the density of both liquids.
The spherical coordinate system in this work is illustrated in Fig.1,with azimuthal angle coordinate φ,polar angle coordinate θ and radial coordinate r.The boundary conditions related to this problem are as follows.

Fig.1.Projected image of the target system and flow direction of continuous phase.
(1)At the droplet interface,the normal velocity is 0:

The tangential velocity is continuous:

The tangential stresses are in balance:

where the shear stress is computed by

(2)At the outer boundary of the field,

(3)At θ =0°and 180°,the velocity vector is continuous[10]:

(4)At φ =0°and 360°,the velocity is continuous:

(5)At the center of the droplet(r=0),the flow is also continuous:

In this study,Eqs.(1)and(2)are solved by a finite volume method in a three-dimensional spherical coordinate system.The computational domain is 0≤r≤R,0≤θ≤πand 0≤φ≤2π,where R is the size of computing domain in the radial direction.At small Re the outer boundary must be larger than the length scale aRe?1/2on which the inertial and viscous terms are comparable.For larger particle Reynolds numbers,the computing domain should cover the area where the vorticity exists.Governing Eqs.(1)and(2)are discretized on a staggered grid,with the nodes for urallocated on the drop surface and the origin.The grid is uniform in azimuthal(φ)and polar(θ)directions,but non-uniform in radial(r)direction.For the internal domain,10-20 nodes are allocated densely and uniformly in the r direction near the surface inside the sphere,whereas away from the surface the nodes are distributed uniformly but with a larger spacing.For the external domain,20-30 nodes in the r direction are set closely and uniformly near the surface since the velocity boundary is very thin,but after that an exponential expansion of cell size is applied:r(n)=r(n ? 1)eα,where α is a small constant used to adjust the node spacing.
The control volume formulation with the SIMPLE algorithm[9]is adopted to solve the governing equations.Although θ =0°and 180°and φ =0°(360°)are boundaries of the computational domain,the fluid flow is continuous there.The values of u and p are specified iteratively as suggested as Zhang et al.[6].Grid sensitivity analysis has proved that R=60a and a grid with[60(internal)+150(external)](r)× 30(θ)× 60(φ)(the minimum Δr=0.0025)suffices for computational accuracy.
The values of the Reynolds numbers would affect the values of velocity and flow structure.The higher the Reynolds number,the larger the shear rate and the velocity near the drop surface.Flow structure will become complex as the value of Reynolds number increases.The streamlines coincide with the trajectories in the steady flow,so the streamlines are obtained by tracing a fluid element.Different topologies of streamlines are obtained around a droplet in the x–y plane with increasing Reynolds number for λ=1.Fig.2 shows that a set of closed streamlines wrap the liquid sphere in the Stokes flow,in which we use Re=10?6as the numerical approximation for Re=0 with zero inertial force.This topology is analogous to that of shear flow field around a solid sphere at Re=0[6].

Fig.2.Streamlines in the x–y plane around a single spherical drop in simple shear flow at Re=10?6 and λ =1.
Fig.3 depicts the topology of streamlines around a liquid sphere at Re=0.1.For such a low Reynolds number(with weak shear), fluid inertia has little effect on the flow field,so that the streamline map is topologically similar to that at Re=0.The trajectory starts from the point(r≈ 1.01,φ ≈ π/2)close to the drop surface.At some instants,a small fluid element leaves the point,spirals outward and goes eventually toward the downstream of the flow field.Thus the streamlines near the liquid sphere are not closed.The smaller the Reynolds number,the longer the spiral pattern.Because of the conservation of mass,a net flux of fluid crosses a concentric circle around the drop in the x–y plane,and this flux is compensated by an inward flux along the vertical direction of the spiral.
As the Reynolds number increases,the velocity becomes larger on both upper and lower sides of the liquid sphere,so does the vorticity,which means that the fluid element here is subject to a force to deform or rotate in the clockwise direction.At moderate Re,the fluid experiences a substantial increase of vorticity along a streamline from the left approaching the top of the drop.The accumulation of vorticity in this region at Re=8 finally results in the topological change of streamlines around the droplet different from that at Re=0.1,as shown in Fig.4.Leal[11]did explain the occurrence of recirculating wakes behind a drop at intermediate Reynolds numbers as a consequence of accumulation of vorticity generated upstream on the body surface.By the symmetry,the same change also occurs on the right hand side of the drop.As in the case of a solid sphere[6],there are wakes with trajectories of reverse direction along the positive or negative x-axis.Some fluid elements from the upstream of flow enter the region near the droplet,proceed to spiral toward the droplet,and eventually leave upwards or downwards away from the x–y plane.As a result of weak vorticity at Re=0.1,there is no wake in the direction of flow.
The flow structure may be influenced by the viscosity ratio,but for different viscosity ratios,the trends of flow field at different Re are similar.
Similarly,the flow field inside a liquid sphere is influenced by the viscosity ratio.The topologies of streamlines in the drop in simple shear flow are obtained for Re=1 and several viscosity ratios.Fig.5 shows that for λ=0.1,the inner liquid is less viscous and the surface shear will intensify the internal flow more drastically.The flow circulations in the droplet are not closed at Re≠0,but the particular topology of streamlines depends considerably on the viscosity ratio.At low viscosity ratios,such as λ=0.1,there are two symmetric points(r≈0.7,φ≈π/2)and(r≈0.7,φ≈ 3π/2),and each is taken objectively as the start-point for all streamlines in the x-y plane.In other words,if any fluid element in the x-y plane moves in the reverse direction,it will reach one of the source points and finally leave the x-y plane.This is a phenomenon for Reynolds numbers from 0.1 to 10.

Fig.3.Streamlines in the x–y plane around a single spherical drop in simple shear flow at Re=0.1 and λ =1.

Fig.4.Streamlines in the x–y plane around a single spherical drop in simple shear flow at Re=8 and λ=1.

Fig.5.Streamlines in x–y plane inside a spherical drop in simple shear flow at Re=1 and λ=0.1.
Fig.6 shows that the topology of streamlines for λ=1 resembles that for λ =0.1 in the x–y plane,but the coordinates of two source points are(r≈ 0.63,φ≈π/2)and(r≈ 0.65,φ≈ 3π/2),which are closer to the origin than those for λ=0.1,and the streamlines are much denser.
At a viscosity ratio of 10,the flow characteristic changes obviously.As shown in Fig.7,there are still two source points on the y-axis,and points(r≈0.3,φ≈π/2)and(r≈0.3,φ≈ 3π/2)are closer to the origin than those for low viscosity ratios.However,the flow element from the source point would reach a closed orbit after spiraling toward the inner surface of the droplet for a certain time and finally leave the x–y plane.The flow field is divided into two parts by this closed orbit.The fluid elements into this orbit starting from the source points circle outwards and eventually reach this orbit,while the outer fluid elements circle inwards and also finally get into the same orbit.In Fig.7,red(inner dumbbell shape)trajectories are the ones passed by two fluid elements from the source points,and the blue(outer)ones by a fluid element from the position(r≈0.8,φ≈0).They merge with each other on the closed orbit.

Fig.6.Streamlines in x–y plane inside a spherical drop in simple shear flow at Re=1 and λ=1.

Fig.7.Streamlines in x–y plane inside a spherical drop in simple shear flow at Re=1 and λ=10.
Although Figs.5-7 represent the changing flow field with several viscosity ratios at Re=1,those flow characteristics are qualitatively consistent with that for other non-zero Re.However,the topology of streamlines,the positions of source points and closed orbits are varied according to the combination of Re and λ.
Subramanian and Koch[7]have proved that the streamlines on the surface for Re=0 are in accord with the Jeffery orbits[13],i.e.,closed curves.This means that every streamline corresponds to the same orbit constant.The Jeffery orbit constant J puts constraint on the coordinates(θ,φ)of the point on the same streamline:

where re=(2/λ+1)1/2is a function of viscosity ratio.When Re is small but not equal to 0,the streamlines at the surface of a droplet are not closed.As a result,the value from Eq.(11)is different for the same streamline but different positions.
Figs.8 and 9 depict the simulated results of streamlines at the surface of a droplet with λ =0.5 for Re=10?6and Re=0.1,respectively.In Fig.8,the streamlines at the surface of a droplet are almost closed(numerical results are suffered from minor errors).For different positions on the same streamline(different θ and φ),the Jeffery orbit constants are very close to the same value.Thus it confirms the analytical results by Subramanian and Koch[12].In Fig.9,one fluid element moves along the trajectory,which starts from a point(θ =0.5,φ =0)at the surface of a droplet.Although the trajectory is densely coiled,we could see that the fluid element spirals toward the x–y plane.It should be noted that the trajectory in Fig.8 is not complete for clarity.A complete streamline at the surface of a droplet should start from the z-axis(θ =0,θ = π/2)and eventually reach the x–y plane(θ = π/2).
The streamlines in Figs.8 and 9 are consistent with the analytical results by Subramanian and Koch[9],who presented only the situation for Re?1.For higher values of Reynolds number,we have to rely on numerical methods.From the numerical results,we can see that the streamlines do not always spiral toward the x–y plane from the z-axis.The pattern of movement relates to Re and λ.
As shown in Fig.10,any fluid element at the surface of the upper hemisphere spirals toward two special points.As the viscosity ratio increases at λ=1,Fig.11 suggests that the fluid element at the surface of a droplet turns to spiral toward a certain closed orbit.In Fig.11,the blue trajectory belongs to one fluid element setting off from the point(θ =0.5,φ =1),while the red one starts from another point(θ=1,φ=1).Finally,these trajectories coincide with each other at the same closed orbit.The fluid element above the orbit runs downward spiral movement,while the one below the orbit spirals upwards.And both eventually circle in the closed orbit.

Fig.8.Streamlines on the surface of a spherical drop in simple shear flow at λ =0.5 and Re=10?6.

Fig.9.Streamlines on the surface of a spherical drop in simple shear flow at λ=0.5 and Re=0.1.

Fig.11.Streamlines on the surface of a spherical drop in simple shear flow at λ=1 and Re=5.
The present numerical analysis reveals more complicated nature of internal,external and surface flows of a neutrally buoyant drop subject to simple shear than that of a drop in simple buoyancy-driven motion.Many drastically different features of fluid flow are demonstrated.These diversified flow structures will certainly influence the mass transfer of a drop in shear flow,and the interphase mass transfer needs further investigation on the basis of resolved flow field.It can be expected that the flow and mass transfer will be challenging for a drop in external flow with shear and convection.
Nomenclature
a drop radius,m
N number of nodes
r dimensionless radial coordinate
R size of computational domain
Re Reynolds number(Re=γ˙a2ρ/μ)
u dimensionless velocity
u dimensionless velocity vector
θ spherical polar angle,(o)
λ interior-to-exterior viscosity ratio
μ dynamic viscosity of continuous fluid,Pa?s
ρ density,kg?m-3
τ shear stress,N?m-2
φ spherical azimuthal angle,(o)
Subscripts
S drop surface
∞ far from the drop
0 initial state or center of a droplet
1 continuous phase
2 drop
Chinese Journal of Chemical Engineering2015年1期