A neutrally buoyant sphere in creeping flow between parallel plates: farfield velocity profiles

A neutrally buoyant sphere in creeping flow between parallel plates: farfield velocity profiles

Chemical Engineering Science, Vol. 45, No. 1, pp. 225-235, Printed in Great Britain. 1990. c ooo%2509/90 S3.W + 0.00 1989 Pergamon Press plc A NEUT...

862KB Sizes 3 Downloads 99 Views

Chemical Engineering Science, Vol. 45, No. 1, pp. 225-235, Printed in Great Britain.

1990. c

ooo%2509/90 S3.W + 0.00 1989 Pergamon Press plc

A NEUTRALLY BUOYANT SPHERE IN CREEPING FLOW BETWEEN PARALLEL PLATES: FARFIELD VELOCITY PROFILES JEFFREY

A. SCHONBERG,

DONALD

A. DREW

and GEORGES BELFORT+ Sciences, Rensselaer Polytechnic Institute, Troy, NY 12180-3590, U.S.A.

Department of Chemical Engineering and Department of Mathematical

(Receised 9 June

1988; accepted

29 March 1989)

Abstract-The farfield creeping flow around a sphere suspended in Poiseuille flow between parallel plates is studied using the method of reflections and the boundary integral method. A sample calculation was performed for the case of a sphcrc located 40% of the channel width from one of the walls. The walls induce crescent-shaped vortices comparatively near but not attached to the sphere, which travel with the sphere. The walls also induce a curious viscous lateral migration totally unlike that observed in Couette flow. In the calculation, the walls diminish the evasive motion of the fluid around the sphere, asymmetrically. That is, the motion on the side of the sphere nearest to a wall is comparatively Battened.

sufficient number to obtain convergence. In contrast, this study treats the sphere as asymptotically small, so that it appears as a singularity to the flow of fluid in the channel. The boundary integral method has also been used to investigate creeping flows. Studies have examined the flow around objects in an infinite medium (Youngren and Acrivos, 1975). This numerical method is based on the viscous hydrodynamic potential theory (Ladyzhenskaya, 1963) which is developed from Green’s theorem. It solves equations to find a surface stress defined on the boundary of the flow instead of solving the Stokes equations directly for functions defined throughout the flow volume. Hence, the number of grid points needed should be less than that required by a finite differencing approach. Furthermore, this method is able to solve problems involving boundaries of general shape whereas Fourier transforms or other “analytical” techniques require highly specialized boundary geometries. On the other hand, the fluid velocities must be calculated in a second step, whereas in a finite differencing approach velocities may be found directly. General information on boundary integral methods can be found in the literature (Belytschko and Hughes, 1983; Fischer and Rosenberger, 1987) In this study, the boundary integral method is “combined” with the method of reflections. The suspended sphere is studied analytically in an infinite medium, an easy task owing to the simple geometry. The walls that bound the channel are analyzed with the boundary integral method. Only one reflection is used but, owing to the shorter-range disturbance associated with a neutrally buoyant sphere, the sphere is not necessarily minuscule. In fact, the expansion, eq. (5) suggests that, for a sphere near the centerplane, the results will be valid to within a few percent for a sphere diameter as large as 10% of the channel width. By removing the sphere surface from the problem with the method of reflections, the boundary integral calcu-

INTRODUCTION

The creeping Bow of a Newtonian fluid is relevant in practical situations with small length scales or high viscosities. When a suspended body enters a flow, the situation is very complex. That body is free to move, unlike the walls of the enclosing channel. Furthermore, the impact of even a small body on the features of a flow such as the pressure drop, though weak, is generally very long ranged. Fluid surrounding a neutrally buoyant sphere does not behave like fluid surrounding a non-neutrally buoyant sphere. Hence, the density of the body is important. The latter situation, which includes gravitational sedimentation, is discussed by Happel and Brenner (1965) among others. Neutrally buoyant spheres have received attention in connection with the lateral migration of spheres in Poiseuille flow (Segre and Silberberg, 1962; Ho and Leal, 1974). Creeping-flow problems including small bodies have been analyzed by the method of reflections (Cox and Brenner, 1967). A sedimenting sphere imposes a long-range disturbance on a flow (proportional to the inverse distance in diameters); however, a neutrally buoyant sphere suspended in a Poiseuille flow spins and imposes a shorter ranged disturbance on the flow (proportional to the inverse distance squared) (Batchelor, 1967). Therefore, the one-term reflection expansion for the latter case is sufficiently accurate for larger spheres than in the former. In this study, the method of reflections is used to analyze the creeping flow around a neutrally buoyant sphere and between parallel plates. Gantos et al. (1980) used a combined analytical numerical

(collocation)

method

to describe the motion

of a sphere and the suspending fluid between parallel plates. The sphere is not assumed to be small. Collocation points are placed on the surface of the sphere in

‘Author to whom correspondence should be addressed. 225

226

JEFFREY

A.

SCHONBERG

lation is simplified, for then the elements of the (wall) surface are planar. The software developed in this study may be used to analyze any flow channel which can be described by planar rectangular surfaces, provided the appropriate version of Poiseuille flow can be calculated. This task is not too difficult in theory if the duct is a straight pipe of arbitrary cross section (Happel and Brenner, 1965). For example, a rectangular duct could be studied. In this analysis, parallel plate geometry is used for its simplicity. The disturbance velocity field is presented. According to the form of the expansion, the results are not valid close to the sphere (within about five radii of the center). However, the fields in this region are easily described since they do not depend on the walls. Hence, the numerical results are intended to describe phenomena which occur comparatively far from the sphere. Integral transform representations of the disturbance field do exist (Ho and Leal, 1974; Vasseur and Cox, 1976) but, to our knowledge, the back transform has not been obtained analytically. MATHEMATICAL

DEVELOPMENT

The origin of the coordinate system is located at the center of the sphere and is moving with the velocity of the ambient Poiseuille flow. The channel Reynolds number, N,,, is defined by the maximum velocity of the Poiseuille flow, U,; the channel width, h (see Fig. 1); and the fluid viscosity. If N,, is much less than 1, the total fluid motion is described by the dimensionless Stokes equations VW conditions v*=V,+nA

For convenience, the disturbance field may be defined as u=v*-u p=p*-P

x=?CIcx*

(3)

where (4) where a is the sphere radius, and h is the channel width. The farfield solution of eq. (1) with eq. (2) for the disturbance without regard to the wall boundary condition is defined as G and is reported elsewhere (Schonberg et al., 1986). It is the Poiseuille analog of the Couette solution (Batchelor, 1967), and may be expressed in the new coordinate as fJi=

-$

xixjx,ejk3c3+O(~5)

where ejk is the rate of strain tensor of the Poiseuiile field (U, P) in terms of the new coordinate

(la)

are xx-

Reflection (ti, @) where

on A

v*=u

on walls

v*-u,

x*-+aZ

is accomplished

velocity,

a new field

ii=u--a

(6a)

j5=p-71n.

(6b)

Note, 7~is the singular pressure field corresponding to m. Hence, the field (ii, p) is regular and satisfies Stokes equations with a non-zero boundary condition at the wall: ii= -e.

and

coordinate

(7)

Since (ii,@) is regular and satisfied the Stokes equations its solution is given by the viscous potential theory (Ladyzhenskaya, 1963). The raw form may be refined by employing the reciprocal theorem over the volume outside the walls. This may be done because both the kernel and the field (r are regular outside the flow channel. Hence ii(x) =

lengths and system.

by defining

(lb)

where the unknown motion of the body is given by the translational velocity V, and its rotational velocity n2,, and where the Poiseuille velocity is U. Note the no-slip condition at the walls is expressed relative to the translating axes, not the laboratory. The body is assumed to be rigid. Since its density is equal to that of the fluid, there is no hydrodynamic force on the body. The body is assumed to have a spatially uniform density. Therefore there is no hydrodynamic torque on it.

Fig. 1. Charkcteristic

(2)

where the Poiseuille flow is given by (U, P). Equation (1) may be easily solved without the no-slip boundary condition at the wall. This boundary condition is recovered by the method of reflections. Specifically, a new coordinate system is defined:

- vp* =o

v.v*=o. The boundary

et al.

c

V(x,q).f(q)

dS.

(8)

This single-layer representation is similar to that given by Rallison and Acrivos (1978). Note that f is not simply the stress at the wall due to (ii, 6) but also

A

neutrally buoyant sphere in creeping flow between parallel plates

involves the stress of the field a. The kernel is Vij(X, q)=

-’

87t

~1

I I[

6,.+(xi-Vi)(xj-Vj)

x-q

‘)

Ix--‘712

1

(9)

.

A corresponding expression for pressure may be written. However, it is not needed to find the velocity field ii(x). Equations (7) and (8) are susceptible to the boundary integral method. The correction term in this expansion-reflection technique is O(K~) so the solution is valid if K2

NUMERICAL

<

1.

TREATMENT

(10)

AND

RESULTS

The numerical solution ofeq. (8) with eq. (7) is based on an approximation of the integral. The bounding surface is discretized and the vector distribution f(q)is assumed piecewise constant. Hence ii(x)=

2 i=l

f(i)

s AW,

X (x,rl)dS

(II)

where the subscript i refers to the ith element of the surface. Once the wall elements are defined, their corresponding integrals are in principle known functions of position x. The unknown “force” distribution may be found by writing eq. (11) for N values of x each restricted to the wall (and placed one per surface element). Equation (7) shows that ii is defined at the wall, given the location of the body. In this way, the “force” distribution f(i) may be found through the solution of an algebraic linear system. Once f(i) is found, the velocity and pressure at any point in the flow domain may be approximated using eq. (11) and its analog for pressure. The accuracy of the approximation hinges on the assumption that the force distribution is piecewise constant. Therefore, in practice, the surface grid must be refined while the velocities are monitored for convergence. The elements of surface AW, are rectangular (of arbitrary size and shape) hence the surface integrals shown in eq. (11) may be evaluated analytically. The analytical solutions however, place some restrictions on the value of x to avoid division by zero. Generally these are well defined geometrically. However, there were unexpected numerical problems on rare occasions. Nevertheless, there were no difficulties in calculating velocities if the surface grid was intelligently designed with a view to the restrictions (Schonberg, 1986). Various difficulties were encountered in the course of the numerical study of the wall-enclosed disturbance field. The first involved the selection of a surface grid to accurately approximate the actual vector distribution f(q). This is not known a priori and so the locations of the areas of high gradients in f(q) are also unknown. The second problem involved memory limitations associated with the Michigan Terminal System operating system which parcels the memory capabilities of the IBM 3081-D computer into small portions. The first of these problems was solved by an

227

educated guess that a grid chosen according to the shape of -a(q), the boundary condition, would be a good first guess. This grid was then refined interactively. The refinement process was as much an art as a science. The logical approach was to subdivide the wall elements where the discontinuity in the vector distribution f(i) was the greatest. This was done but another factor was important. Larger wall elements had a much greater influence on the system than did smaller elements. Thus after a point the larger elements had to be subdivided even though the discontinuities with their neighboring elements were not the largest in the system. In the course of the refinement process the outer boundaries of the grid were set by monitoring the influence of the outermost surface elements on the velocity field u. No attempt was made to take advantage of asymptotic forms of the function f. The refinement of the grid could have been more systematic if not for the memory limitations which restricted the number of wall elements. The operating system effectively limited this number to 100. The results presented were developed from 360 wall elements by taking advantage of symmetry in the parallel plate boundaries and in the fundamental solution Vij [eq. (9)]. The following relationships exist: If f=(f,,fi,S3)

then f’=(-fi,fi,fs)

and

where f is the “force” on the element, f’ is the “force” on the element which is a mirror image through the 2-3 plane, and f”, is the “force” on the element which is a mirror image through the l-3 plane. The velocities u satisfy the same relationships. Thus the number of wall elements required was reduced by a factor of 4. The problem of storing the matrices could have been avoided if the linear system was solved iteratively. Gauss-Seidel iteration was investigated but convergence is only guaranteed if the matrix is diagonally dominant. In another application of the boundary integral method (Youngren and Acrivos, 1975) the matrix is reported to be not diagonally dominant. This would appear to be true for this system. The crucial factor is the long range of the singular solution vij The finer the grid is made the worse the situation becomes. As previously mentioned, convergence was determined by comparing velocities at a large number of locations in the flow domain for several grids in a series (the final iteration was checked at 405 locatations). The densities f were not converged because of the memory restrictions. A broad sample of velocities was taken because grid refinements often had a regional, but not a global impact. The velocities (which range in magnitude from zero to 40) had settled to within 0.03 with good confidence. The fields are clearly order-one quantities as is shown in the graphical results presented in this section. Hence the convergence is good. The results are order one because the K’

228

JEFFKEY

A.

SCHONBERG

factor in cr may be taken out of the linear calculation. Hence, the numbers given must be multiplied by ~~ to give the velocity scaled by U,, the maximum Poiseuille velocity. Although convergence is good for these velocities, meeting the “no-slip” condition between the “zero” points is much more difficult. Hence, the program was not used to study velocities extremely close to the wall. The results are shown as isovelocity line contour maps in Figs 24. In each case the lower wall is located at (0, O,-0.4). These show the disturbance velocity fields in the first quadrant (with respect to the l-2 axes). The Poiseuille flow is in the 2-direction and varies in the 3-direction according to U=e,(0.8x3-4x5).

(12)

et al.

The disturbance velocity fields in the other three quadrants are constructed through the symmetry rules. Due to anti-symmetry (and continuity) certain zero-contour lines exist by definition. In Fig. 2, eight maps of u1(x,, x3) are shown for eight values of x1. In each of these the x,-axis is a line of zero velocity. (The x,-axis and the line x3 =0.6 are also zero-velocity contours by virtue of the no-slip condition.) Had a map been drawn for the case of x, equal to zero, velocity would have been zero everywhere. In Fig. 4, eight maps of u,(xl,xZ) are shnwn for eight values of xj. In each of these the x,-axis is a zero-contour line. Incidentally, in Fig. 3 the x,-axis and the line xX =0.6 are zero-velocity contour lines due to the no-slip condition. Note that, in Fig. 4, the sphere is on the x,-axis.

Ul(X2,

06.

Fig. 2.

(a).

x3)

x, - ,325

A neutrally buoyant sphere in creeping flow between parallel plates

UI(XZ, x3)

Ul(X2,

x3)

x, = -425

UI(X2, x3)

x, = .8

Ul (x2. x3)

229 x, = .65

x, = 1.0

06

06

Fig. 2. (a) Isovelocity contour lines for disturbance flow in the l-direction, for x1 =0.04, 0.125, 0.225 and 0.325. Negative and positive values represent flow into and out of the page, respectively. Cross hairs are not locations of extreme velocities. (b) Isovelocity contour lines for disturbance flow in the l-direction, for x, = 0.425,0.65,0.8 and 1.0. Negative and positive values represent flow into and out of the page, respectively.

DI!XXJSSION Overview

of results

The results are useful for studying particle interactions and for studying the disturbance velocity itself. Although the disturbance velocity is only a small portion of the total velocity field, it may be considered as a separate entity since the creeping-flow equations are linear. The wall-free disturbance field consists of currents which flow radially inward or outward from the sphere. These currents are driven by a mass conservation effect. The fluid surrounding the rigid rotating sphere rotates and strains, changing its shape. These currents do not exchange fluid between them, as

shown by eqs (5) and (12). Specifically, the three planes defined by the coordinate axes experience no normal velocity at any location. Thus no fluid passes between the octants. Global continuity is achieved through the decay of these current at infinity. The presence of walls will obviously alter this feature of the disturbance flow. However, some features are the same. Figures 2Z4 show strong currents close to the sphere (on the channel scale) which are roughly symmetric in magnitude and are directed according to eqs (5) and (12). It is not surprising that the wall effect is fairly weak near the sphere where the rigid body

effect is strongest.

On the other

hand

the

230

JEFFREY A. SCHONBERG

contour graphs clearly show the influence of the walls. Fluid passess across the 1-2 plane (Fig. 4), and across the l-3 plane (Fig. 3). Fluid does not pass through the 2-3 plane due to the symmetry in the flow channel in the l-direction. Hence the walls cause the outflows and the inflows to feed each other, as shown in Fig. 5. Although the disturbance flow has closed loops like separation flows, the total flow does not because the disturbance flow accounts for only small direction changes in the Poiseuille flow. However the disturbance velocity contours shown in Fig. 4 (x3 = 0.0) occur in the plane of the particle where the undisturbed velocity is zero. Hence the updraft which occurs in the region downstream of the sphere corresponds to the redirecting of some of the flow from the oncoming Poiseuillean current below the sphere, in the xj sense,

x3)

U2(Xlt

et al.

to the overtaking Poiseuillean current above the sphere (see Fig. 6). [Note a similar effect occurs in the plane x3 =0.2 (Fig. 4) where the undisturbed velocity is also zero. The result is an eddy. Due to the velocity in the l-direction (Fig. 2) this exchanges fluid with the rest of the flow and is not a separated flow.] Corresponding to the updraft a mirror image effect occurs upstream of the sphere also. The final feature to be pointed out in the disturbance flow is the damping of the third component of the velocity by the bottom wall of the duct. Ifthere was no wall effect, the isovelocity contours in Fig. 4, in the xx = -0.1 and x3 = 0.1 planes, would be identical. Likewise the xs = -0.2,0.2 maps would be identical. However, the velocities in the 0.1 and 0.2 planes are noticeably larger.

x2 = .a4

x2 = ,125

U2(XI> x3)

_I 0.2 (a)

04 -2.0

-0.5

-0.2

0.0

x3

-02 -

00

x3)

UZ(Xi.

x2 = .225

uz(x7,

__

04

‘-K

-0.2

-1 .o

-5.0

02

\\

-10.7

x3

t

00

02

_J

+

Fig. 3. (a)

04

x3)

06

08

xg = .325

10

A neutrally buoyant sphere in creeping flow between parallel plates U2(Xll

x3)

x2 = .425

UZ(XI.

231 xz = .65Q

x3)

06 (b) 04 --

-0.0

--

-02 -

02

I 06

04

I 08

r 10

x1

Fig. 3. (a) Isovelocity contour lines for disturbance flow in the 2-direction, for xz =0.04,0.125, 0.225 and 0.325. Negative values represent flow out of the page. Cross hairs are not locations of extreme velocities. (b) Isovelocity contour lines for disturbance flow in the 2-direction, for xI =0.425, 0.650.0.8 and 1.0. Negative

values represent flow out of the page.

Viscous lateral m&p-ations The viscous lateral migrations (x,-direction) of two or more spheres is easily determined from the singlesphere disturbance field. On the scale of the channel, the spheres appear as points to first approximation. As the spheres experience no body forces they move with the fluid. Furthermore, since they are represented as points, they have no second-order interactions. The first-order interactions are found directly from the contour maps. The symmetry results discussed above imply that the lateral velocities (u3) downstream of the sphere are the reverse of those upstream. Since this is true for spheres at any location, there wiIl be no net

lateral migration of a sphere which experiences a flow ofaxially uniform concentration. However, a sphere in the leading edge of a cloud of spheres will migrate. This result is qualitatively the same as the wall-free results (Schonberg et al., 1986). The analogy breaks down in the details. Figure 4 shows that two possible lateral interactions occur due to the wall effect, whereas in the wall-free flow only one interaction may occur. These are quantitatively diagrammed in Fig. 7. Figure 4 shows that the wall effect may induce migrations in situations where none would otherwise occur. In general, the lateral velocities are strongest near the particle (on the scale of the channel) where they are

232

er al.

JEFFREY A.SCHONBERG

most

like the wall free flow. The “reverse-flow”

velocities velocities.

lateral

are always about 10% of the “core flow” The qualitative nature of axial dispersion of

well-separated particles should be the same in the general case as in the wall-free case. The main difference is that the migrations in the general case are somewhat weaker. Finally, note that the particles have an inertial migration velocity of order N,,K~, where N,, is the channel Reynolds number defined previously (Cox and Brenner, 1968).

U,(Xl.

x2)

This velocity is therefore much smaller than that discussed here but may be important in situations in which the symmetries in the viscous fields cause their effects to cancel. Vortices The vortices predicted here are probably a device to reduce viscous stress in the fluid. In a sense, they may act as roller bearings, easing the motion of the fluid around the sphere (see Fig. 8). These vortices are quite narrow, their width being of magnitude lc3.

XJ =-a3

U,(Xl>

x2 = -0.2

x2)

IO-

10

(4 08

08

-

-\\\

XP

O6 0.2

0.1

x2

Ir

Qy, 00

0.2

0.5

06

,,,, I

02

06

04

06

10

00

x1

02

04

06 x1

Fig. 4. (a).

08

10

A neutrally

U3(Xl. 1

x2)

buoyant

sphere in creeping

flow between parallel

xg = 0.4

plates

233

x3

=

0.5

n-

“2

Fig. 4. (a) Isovelocity contour lines for disturbance flow in the 3-direction, for xj = -0.3, -0.2, -0.1 and 0.0. Negative values indicate flow into the page. Cross hairs are not locations of extreme velocities. (b) Isovelocity contour lines for disturbance flow m the 3-direction, for xj =O.l, 0.2, 0.4 and 0.5. Negative values indicate Row into the page.

CONCLUSION

This work demonstrates the use of the method of reflections in concert with the boundary integral method for a sample creeping-flow problem. By taking advantage of the inherent symmetry of the system, the memory requirements were reduced by a factor of 16. This allowed direct solution of the linear algebraic system. The results may be used to study first-order interactions between neutrally buoyant spheres. The size restriction on such spheres is less extreme than for sedimenting spheres. In principle, more general prob-

lems involving non-spherical particles and more complicated ducts could be treated by this algorithm. A major improvement would be to avoid the need to solve the system of linear algebraic equations directly, which requires the storage of large dense matrices.

Acknowledgements-The authors wish to acknowledge the support of the National Science Foundation for funding this work under grant number CEP83-14443. One author (DAD) also acknowledges the support of the U.S. Army Research Office under contract DAAG 29-85-G-0088.

234

JEFFREY A. SCHONBERG et al.

b)

b) Fig. 7. Two kinds of evasive maneuvers possible: (a) evasive maneuver similar to that found in infinite Couette and infinite Poiseuille flow (wal free). (b) Evasive maneuver due to wall effect (the reference frame in each case is between the particles, travelling with Poiseuille flow).

Fig. 5. (a) Qualitative description of disturbance flow relative to the sphere, side view (Poiseuille flow is from left to right and defines the term “downstream”). (b) Qualitative description of only the largest downstream loop of the disturbance flow, relative to the sphere, shown in (a) (viewed from downstream of the sphere).

Fig. 6. Qualitative description of the total creeping flow, relative to the sphere with exaggeration of lateral velocities. The four recirculation flows are actually narrower than the sphere, their width being ICY,if the channel width is 1.

NOTATION ejk

e2

f

f(i) h N,,

P

rate flow

of strain tension of (dimensionless) given

the Poiseuille by eq. (5)

Fig. 8. Detached vortices. The vortex shape is determined the 0.0 contour in Fig. S(b) (x, =0.2).

P*

fluid

P

Poiseuille pressure (dimensionless) dimensionless disturbance field

u

pressure

by

(dimensionless) defined

unit vector in the 2-direction (see Fig. 1) unknown vector function in integral eq.

by eq. (2) Poiseuille

(8) discretized version of f refering to the element of surface channel width channel Reynolds number based on the maximum velocity of the Poiseuille flow and the channel width dimensionless disturbance pressure defined by eq. (2)

maximum velocity of the Poiseuille flow fluid velocity (dimensionless) Kernel in the integral eq. (8) corresponding to the creeping flow ,equations given by eq. (9) translational velocity of sphere A (dimensionless)

V* W

AWi

pressure

(dimensionless)

surface defined by the duct ith element of surface

walls

A neutrally x

magnitude

X

dimensionless

X*

dimensionless

xi,

the sphere x expressed

Greek

xk

sphere

in creeping

vector,

scaled

by

vector,

scaled

by

of x

the channel

xj>

buoyant

position width position radius in tensor

notation

letters position vector variable of integration (dimensionless) q expressed in tensor notation ratio of sphere radius to channel width pressure singular

field corresponding to solution of dimensionless

ing flow equations CTexpressed rotational

creep-

given by eq. (5)

in tensor notation velocity of sphere

A (dimen-

sionless) Laplacian

operator

gradient operator vector dot product much ordering

less

f(x)=Cs(x)l

parallel

plates

defined

if lim

x-x o

f(x) g0

by is finite

REFERENCES

Batchelor, G. K., 1967, An Introduction to Fluid Dynamics Cambridge University Press, Cambridge. Belytschko, T. and Hughes, T. J. R. (Eds), 1983, Computational Methods for Transient Analysis. NorthHolland, Amsterdam.

235

Cox, R. G. and Brenner, H., 1967, Effect of finite boundaries on the Stokes resistance of an arbitrary particle, Part 3, translation and rotation. J. Fluid Mech. 28, 391411. Cox, R. G. and Brenner, H., 1968, The lateral migration of solid particles in Poiseuille flow--l. Theory, Chem. Engng Sci. 23, 147-l 73. Fischer, T. M. and Rosenberger, R., 1987, A boundary integral method for the numerical computation of the forces exterted on a sphere in viscous incompressible flows near a plane wall. J. appl. Math. Phys. 38, 339-365. Gantos. P. R.. Pfeffer. R. and Weinbaum. S.. 1980. A strong interaction theory for the creeping motion of a sphere between plane parallel boundaries, Part 2, parallel motion. J. Fluid Msch. 99, 755-783. Hap@, J. and Brenner, H., 1965, Low Reynolds Number Hydrodynamics, pp. 33-35. Prentice-Hall, Englewood CM, NJ. Ho, B. P. and Leal, L. G., 1974, Inertial migration of rigid spheres in two-dimensional unidirectional flows. J. Fluid Mech. 65, 365400. Ladyzhenskaya, 0. A., 1963, The Mathematical Theory of Viscous Incompressible Flow. Gordon & Breach, New York. Rallison, J. M. and Acrivos, A., 1978 A numerical study of the deformation and burst of a viscous drop in an extensional flow. J. Fluid Mech. 89, 191-200. Schonberg, J. A., 1986, Thesis in partial fulfillment Rensselaer Polytechnic Institute, Troy, NY.

than operator

flow between

of PhD.

Schonberg, J. A., Drew, D. A. and Belfort, G., 1986, Viscous interactions of many neutrally buoyant spheres in Poiseuille flow. J. Fluid Mech. 167, 415426. Segre, G. and Silberberg, A., 1962, Behavior of macroscopic rigid spheres in Poiseuille flow, Part 2: experimental results and interpretation. J. Fluid Mech. 14, 116-134. Vasseur, P. and Cox, R. G., 1976. The lateral migration of a spherical particle in two-dimensional shear flows. J. Fluid Mech. 78, 385413. Youngren, G. K. and Acrivos, A., 1975, Stokes flow past a particle of arbitrary shape: a numerical method of solution. J. Fluid Mech. 69, 377403.