Slow Variations in Continuum Mechanics

ADVANCES I N A P P L I E D M E C H A N I C S , V O L U M E

25

Slow Variations in Continuum Mechianics MILTON VAN DYKE Departments of Mechanical Engineering and Aeronautics & Astronautics Stanford University Stanford, California 9430.5

I. Introduction Any good plumber can tell us that the velocities across a 1-in. pipe are just four times those in a 2-in. pipe to which it i s connected (or at any rate they would be if those nominal diameters were the actual ones). He might concede, however, that this common sense concliision becomes questionable in the vicinity of the juncture of the two pipes. This hydraulic approximation is a simple and useful idea which, as the quasi-cylindrical approximation, has its counterpart in other branches of mechanics. It has long been applied to specific problems in many research papers and textbooks, but no general exposition of it can be found. Here we undertake to discuss it from a unified point of view, with illustrations drawn mostly from fluid mechanics and the linear theory of elasticity. We recognize that the approximation is founded on the assumption that the boundaries of the region being considered vary much more slowly in some directions than in others. We call this a slow variation. (We could equally well regard it as a rapid variation in the transverse direction, because what matters is only the relative rate of variation of the geometry.) In three dimensions we can distinguish two classes of slowly varying geometry. The more common type is that of a .rfender object, exemplified by a needle, whose two transverse dimensions vary slowly in the longitudinal direction. The other type is that of a thin object, such as a razor blade, whose thickness varies slowly in the other two dimensions. Even in a slowly varying domain the solution could vary rapidly, as for waves traveling along a slowly varying channel. We exclude this important 1 Copynght 0 1987 by Academic Press, Inc All rights of reproduction m any form reserved

2

Milton Van Dyke

situation, to which a great deal of interesting work has been devoted. Thus we consider problems of equilibrium, governed typically by elliptic partial differential equations such as Laplace’s equation. (Note that from this point of view, steady fluid flow is an “equilibrium” situation.) We also disregard ends, edges, junctures, and other discontinuities. These exert only a local influence, which in most cases decays exponentially and so becomes negligible beyond a few widths. Viscous flow is an exception: although a disturbance can propagate only a few stream widths upstream, it is swept downstream and so extends in that direction a distance of the order of the width times the Reynolds number (Van Dyke, 1970). Local solutions for such discontinuities can be joined to the slowly varying solutions that we shall consider, using the method of matched asymptotic expansions, to render them uniformly valid. Thus Keller and Geer (1973, 1979) have connected thin slowly varying streams of a heavy inviscid fluid in a free jet, along a solid wall, or in a channel. We seek to embed the intuitive idea of slow variations into a systematic scheme of successive approximations. The result is a regular perturbation expansion of the solution in powers of a parameter that characterizes the slow variation. We can in principle extend the series to arbitrarily high order. In practice, of course, the computational labor grows so rapidly that only a few terms can ordinarily be calculated by hand; but it has increasingly been found possible to delegate that labor to a computer, using either purely arithmetic programs or the newer ones that carry out symbolic manipulation. The question then arises whether the series actually converges for at least some range of the perturbation parameter or is merely asymptotic with zero radius of convergence. We endeavor to answer this question both analytically and by study of numerical results. The key to treating a slow variation is to rescale the coordinates in different directions so that the variation formally becomes “normal.” This transfers the perturbation parameter from the boundary conditions to the differential equations, which can accordingly be simplified by approximation. It is remarkable that this simplification also renders the approximate solution uniformly valid. That is, without the rescaling, the perturbation solution is in many cases valid only for slight variations and breaks down where the variation is appreciable. With rescaling, on the other hand, arbitrarily large variations are acceptable, provided that they are slow. In this way our plumber can accurately treat a 1-in. pipe that expands to a 10-in. pipe, provided that it does so slowly. Thus a slow variation may in its original form be regarded as a singular perturbation problem, and the rescaling as the simplest technique for transforming a singular to a regular perturbation. This procedure was perhaps first systematically applied by Blasius (1910) to the steady plane laminar flow through a slowly varying symmetric channel. It is remarkable that within fluid mechanics it has subsequently been extended to other viscous flows, whereas the simpler potential flows have

Slow Variations in Continuum Mechanics

3

been largely neglected. We take advantage of this gap by first solving Laplace’s equation in many of our illustrative examples. 11. Two-Dimensional Shapes

The method of slow variations has been applied mostly to thin twodimensional shapes in the plane or slender axisymmetric ones in three dimensions. The governing partial differential equations can then be approximated by a succession of ordinary differential equations, and in fact in most applications the problem is reduced to mere quadratures. We consider such examples here, deferring fully three-dimensional problems to subsequent sections.

A. SYMMETRIC PLANESTRIP Consider a n infinite symmetric plane strip of slowly varying width, such as that sketched in Fig. 1. By slow variation we mean that although the width of the strip may change considerably, it does so slowly, because the slope of the boundary is small, say, of order E << 1. This may be expressed formally by writing the equation of the boundary as y = *f(Ex). Here x and y have been made dimensionless by reference to some characteristic transverse dimension of the strip, such as its minimum or average width. We assume that the derivative f’is of order unity, as are the combinations ff”, f 2 , f”’, and so on. 1. Potential in Hyperbolic Strip

To illustrate the essential simplification associated with slow variation, we consider the thin symmetric strip that is bounded by the hyperbolas y = *(l + E ~ x * ) ” * . In this example the width actually tends to infinity in both directions, but it does so slowly. Suppose that we seek within the strip

FIG. 1 .

Slowly varying symmetric strip.

Milton Van Dyke

4

a solution 9 of the Laplace equation lcIXx + $ ,, = 0 that assumes the values *l on the upper and lower edges. Here could be the steady temperature distribution in the strip with its edges maintained at different constant temperatures or the stress function for irrotational flow through a hyperbolic channel and so on. The exact solution of this problem can be found in elementary closed form using elliptic coordinates (Morse and Feshbach, 1953, pp. 1195-1196). Instead, we seek an approximate solution when E is small. We carry out a perturbation solution in powers of E in two different ways. First, we ignore the slowness of the variation of width and calculate a conventional perturbation expansion leaving the longitudinal and transverse coordinates unaltered. Second, we recognize the slow variation by rescaling those coordinates relative to each other before introducing the perturbation expansion. It is illuminating to compare these two different approximations, which we characterize as slight versus slow variations. a. Slight Variations In the conventional first method, in terms of the original coordinates x and y, the strip approaches a strip of constant width 2 as E tends to zero. For that the solution is simply (I, = y. Then expanding the equation of the boundaries in Taylor series for small E yields only even powers, which suggests that the perturbation solution will also proceed in powers of E ' , in the form (I,(% y ;

=Y

+ E 2 9 2 ( X , y ) + E493(X, Y ) + . . . .

(2.1)

Substituting this expansion into the Laplace equation and equating like powers of E' shows that each term separately must satisfy that equation. It is a familiar feature of perturbation methods (Van Dyke, 1975, Section 3.8) that when, as here, the position of a boundary varies with the perturbation quantity E , it is necessary, in order to carry out a systematic expansion scheme, to transfer the boundary condition by series expansion to the basic position of the boundary corresponding to E = 0. Here that transfer can be effected by Taylor series expansion. Thus the requirement that (I, = 1 on the upper edge of the strip is expressed in terms of the successive (I,, and their derivatives evaluated at y = 1, and similarly for the lower edge. Then like powers of E~ can be equated to give the boundary conditions and so on. Plane harmonic functions satisfying these conditions are readily found as polynomials in x and y (or more easily as the real parts of odd polynomials in the complex variable x + iy). The resulting perturbation series is (I, = y [ 1 - i E Z ( 1

x (43

+2

+ 3x2

-

y2) +

A&"

+

1 0~ 70y2 ~ + 1 3 5 -~ 2~ 7 0 ~ ~ 27y4) ~ '

+

*

. .].

(2.3)

Slow Variations in Continuum Mechanics

5

If this is interpreted as the stream function for inviscid flow in a channel, the derivatives of CC, give the velocity components. Thus the axial velocity is U =

$!ly

+ 3X2 - 3y’) + + 2 1 0 ~ ’- 2 1 0 +~ 1~3 . 5 -~ 810x2y2 ~ + 1 3 5 ~+~ .) . . .

= 1 - ;&*(I

x (43

j&E4

(2.4)

This approximation evidently fails for large x. Whether we examine CC, itself or its derivatives, the second term becomes as large as the first when x is of order 1 / ~ and , that difficulty is compounded in the third and higher approximations. It is clear physically why the approximation breaks down at such great distances. When x is of order 1 / ~the , boundary has departed significantly from the basic uniform strip, and the boundary condition has to be transferred over too great a distance to remain valid. b. Slow Variations

The preceding procedure is not the usual practical way of calculating the flow in a channel. Our intelligent plumber would instead apply the hydraulic approximation, assuming that at any station the velocity is parallel to the axis and constant across the channel. This gives the familiar quasi-onedimensional result that the speed u varies inversely as the cross-sectional area, in order to satisfy the requirement of continuity. In our special case of the hyperbolic channel the result is

(2.5) This hydraulic approximation has the great advantage over our preceding result that it is not restricted to slight variations of channel width. It is accurate even for enormous variations, provided that they take place slowly enough. Furthermore, it soon recovers its accuracy outside any local zone of rapid variation, such as a juncture or discontinuity. This is a special simple example of the quasi-cylindrical approximation, which is familiar in many other problems and fields. Some examples are current flow in a slowly varying wire and bending or torsion of a slowly varying beam or shaft. On its face the quasi-cylindrical approximation is an ad hoc rather than a systematic one. That is, it is not immediately obvious how it can be embedded into a rational scheme for successively calculating higher approximations. In fact, however, this is easily accomplished by first shifting the perturbation from the boundary conditions to the differential equation. This is accomplished simply by stretching the coordinates differentially; that is, we “square up” the geometry by introducing the contracted abscissa X = E X . Then the small parameter e disappears from the equation of the boundaries, which becomes y = f ( 1 + X 2 ) ’ ” . [Only the relatiue stretching of longitudinal and transverse coordinates is significant. Hence instead of contracting the longitudinal coordinate, we could alternatively magnify the transverse one by introducing y / e in place of y , and we would prefer this

6

Milton Van Dyke

choice when the original problem is scaled such that the region is “thin” rather than “long”-for example, if we had defined the hyperbolic strip on a smaller scale according to y = f ~ ( 1+ x2)’/’.] With this distortion of the coordinates, the perturbation parameter appears in the differential equation rather than the boundary conditions, as the problem becomes* E2+hxx

+ Gyy = 0,

*

at

= Zt 1

y

=

*(l

+X2)1/2.

(2.6)

As E tends to zero the boundaries now remain hyperbolic rather than degenerating to straight lines as in the preceding approximation. On the other hand, the differential equation does not remain Laplace’s equation, but reduces to simply t,byy = 0. Hence elementary quadratures yield the first approximation $1

=-

Y

G-2’

In the application to inviscid flow through a channel, this stream function gives the longitudinal velocity (2.5) of the hydraulic approximation. Furthermore, working with the stream function in this way has the advantage, because continuity has been satisfied exactly, that it also provides a first approximation to the transverse velocity, which, being small, is absent from the hydraulic approximation

(Here in the last form we have restored the original abscissa, the contracted one having served its purpose.) Had we treated the flow problem using the velocity potential rather than the stream function, this result for the transverse velocity would have been delayed to the second approximation. We can now systematically calculate higher approximations by assuming a regular perturbation series. Again the appearance of .s2 in the problem suggests an expansion in even powers, as

*=Jr +&2*,(x,y)+s4$3(x,y)+.... + x2

(2.9)

Substituting into (2.6) gives the problem for the nth approximation as Thus at every stage, instead of seeking a solution of Laplace’s equation, as in the previous approximation of slight variations, we need merely perform * A mathematician would also change the notation for $ to indicate that it is not the same function of X and y as it was of x and y, but we follow the usual engineering practice of using the same symbol when no confusion can result.

7

Slow Variations in Continuum Mechanics two quadratures. Thus we easily find the second approximation as

There is now no symptom of nonuniformity in this or any higher approximation; the perturbation is regular. For example, far downstream our hyperbolic channel approaches asymptotically the wedge-shaped channel y = + EX = + X . There our second approximation (2.1 1) becomes s=v[,+iE2(l

-$)+

1

X

0(E4,x-2)],

(2.12)

and this is just the expansion to order F * of the exact solution for a wedge of semivertex angle tan-' E, which may be regarded as produced by a source at the vertex:

*

=

tan-'( & y / X ) tan-'(y/x) tan-' E tan-' E

(2.13)

Expanding this in Taylor series for small E reproduces our approximation of slow variations (2.12) and shows that it converges for E < 1, that is, up to a wedge half-angle of 45". The contracted abscissa X having served its purpose, we can now restore the original coordinate x in o u r approximation (2.11) for the hyperbolic channel to obtain

*=d- 1 +

1 - 2E'X' y2 E2X2 [ l - E ' 6(1 E ~ x ' ) ( '-IfE2x?)

+

]

+ O(c4) .

(2.14)

As indicated by the final order symbol, our result in this form can still be regarded as an expansion in powers of E'. Now, however, the coefficients also depend explicitly on E. Hence this is no longer a power series or even an asymptotic expansion in the classical sense. ErdClyi (1961) calls such a series a generalized asymptotic expansion. If we further expand this result formally for small E, we reproduce our previous result (2.3) of slight variations. This process makes clear how the validity of the solution is thereby destroyed for large x. More precisely, it shows that whereas our result (2.14) of slow variations converges for E < 1, that of slight variations converges only for EX < 1 . It may seem remarkable that the approximation of slow variations is not only much simpler to calculate than that of slight variations-involving only quadratures rather than integration of partial differential equationsbut at the same time yields more complete results. It must be understood, however, that these twin advantages arise only when the variation is actually slow. When the shape varies on a normal scale rather than slowly (Fig. 2), only the first method can be applied.

8

Milton Van Dyke

FIG. 2.

Slightly varying symmetric strip.

If we interpret t,b as the stream function for potential flow, the variation of speed along the axis of the hyperbolic channel is as shown in Fig. 3, with E = 0.8. Even for this rather large value of the perturbation parameter, the first approximation of slow variation (the hydraulic approximation) yields reasonable accuracy throughout, and the second approximation lies considerably closer to the exact result. The latter is found, using elliptic coordinates, as u(x, 0) = &(tan-' &)-'[I

+~

~

+ (x ~1) ] - ' ' ~ .

(2.15)

This shows that, as for the limiting wedge, the approximation of slow variation converges for E < 1. 2 . Potential in General Symmetric Strip The approximation of slow variations is so simple for a symmetric strip that it is actually easier to treat the general shape than any particular one.

FIG.3.

Speed on axis of hyperbolic channel in potential flow.

Slow Variations in Continuum Mechanics

9

Thus for the strip of Fig. I , described by y = f f ( ~ x ) we , easily find the third approximation of our potential problem as

cc, = (Y/f)+ &'(l/f)"Y(f'

Y') + &4{&[f2(1/f)"1"Y(f2- y 2 ) - M / ~ ) ( ~ w ~- y4)) + o(E6). (2.16) -

This expansion terminates at the first term for the unrealistic shape with f ( X ) = (A + BX)-' but is otherwise an infinite series. In the corresponding viscous problem discussed later, in Section 4, Lucas (1972) caries out the indicated differentiation and replaces the transverse coordinate y by the fractional distance 77 = y / f across the strip. Here that gives the more explicit but lengthier form

cc, = 77 +;&'(2f2-ff)(7/

+ 4f'f'f"' + 8f2f'p

+

- f3f'4')(

lOff'S"+ 77 - v3)- &(24ff4 - 36ff"f"

- f3f4')(

7 - T ~ ) ] O( E

-

773)

&4[&(4f'4-

+

4f2p

+ 6f2Y2 (2.17)

~ ) .

This illustrates that the term in E " (where n is even) involves the first n derivatives of the shape function f: More precisely, it is a polynomial in odd powers of 77 of degree n + 1 multiplied, as Lucas (1972) observes for the corresponding laminar-flow solution, by a linear combination of all possible n-tuple products o f f and its first n derivatives f, f',f", . . . ,f'"' that involve n differentiations. 3 . Elastic Strip in Tension We consider now a slightly more complicated problem, from the linear theory of elasticity. Let the symmetric strip of Fig. 1 be subjected to tension at its distant ends. It is in a state of plane stress if it is very thin and plane strain if it is very thick. The problem for the stresses is the same in either case. It is convenient to work with Airy's stress function 4, according to which the three components of stress are given, in the notation of Timoshenko and Goodier (1951), by ax

=

4,,,

fly

=

A x ,

Txy

= -4xy.

(2.18)

Then 4 satisfies the biharmonic equation. In terms of the contracted abscissa X = E X , this becomes

,+,

2

+ 2.5 4xx,, + E 4 4 X X X X

=

0.

(2.19)

We see that again the approximation of slow variations reduces the partial differential equation to a succession of quadratures. The conditions that the edges of the strip y = * f ( X ) are free of stress become

Milton Van Dyke

10

If all variables are made dimensionless by reference to the tension force and a characteristic transverse length, the condition that the prescribed tension acts across each section of the strip is

4y

= i2

at Y = f ( X ) ,

(2.21)

and 4 is an even function of y by symmetry (so that there is no bending moment in the strip). For E = 0 the differential equation degenerates to ~ y v y 4 .= 0, and then quadratures yield the first approximation for the stress function

4

= t[Y2/f(X) + A X ) ] .

(2.22)

Hence the stress components are

The first of these is the simple quasi-cylindrical approximation according to which the tension is uniform across each section; but again the use of an auxiliary function has yielded a bonus, giving first approximations also for the other two smaller components of stress. Seeking corrections of order E’ and E~ to the stress function yields, after some calculation, the third approximation

4

= t(Y’/f+f)

+

kE2(l/f)”(f2 - Y2I2 - -ii[fz(i/f)”]pt + &,(1/f)(4)(2f2

-

E4{-&f(4)

+ y2))(f2 - y2)2+ q E 6 ) . (2.24)

This has a structure resembling that of the potential solution (2.16) and reflects the kinship of the harmonic and biharmonic equations. That resemblance is enhanced by writing this result in terms of the fractional distance 7 = y/f across the strip (and carrying out the differentiations), which gives

4

= tf{(l

+ 7’)

- iE2(2f2 -ff”)(l

-

2T2 + q4)

+ - =&(8f4 - 20ff”f” + 8f2y2 + 8 f 2 f , y + f’f‘“’) + A(24f” - 36ff”f” + 6f2f,r2 + 8f’fP - f’f‘”’, E ~ [

x (2 + q2)](1 - 72)2)+ 0 ( E 6 ) .

(2.25)

We find it informative to compare the solutions of our various examples for similarities of structure in their dependence on the perturbation parameter E, on the widthf of the strip and its derivatives, and on the fractional ordinate 7.Here, comparison with the corresponding form (2.17) for Laplace’s equation shows that, aside from a factor of 7 there and o f f here, the series solutions of the harmonic and the biharmonic equation have the same structure in that they both proceed in powers of e2, and the nth terms

Slow Variations in Continuum Mechanics

11

contain the same combinations of the ordinate f and its derivatives, multiplied by polynomials in v 2 that are of degree n for the harmonic equation but degree ( n + 1) for the biharmonic equation. The special case of a wedge (f = X ) can be solved in simple closed form (Timoshenko and Goodier, 1952, Section 35). In polar coordinates the stress function is

4

=

( r e sin 0)/(2a

+ sin 2 a ) ,

(2.26)

where the semivertex angle is a = tan-’ E. Expanding this in powers of E checks the result from our third-order slowly varying approximation (2.24)* and shows that again it converges for E < 1. 4. Laminar Flow in Channel

The approximation of slow variations has been carried farthest in the problem of steady laminar flow through a symmetric plane channel. Two different schemes of successive approximation have been proposed, depending on the magnitude of the Reynolds number R. If the Reynolds number is of order unity, each approximation requires only quadrature of polynomials, as in the preceding problems. However, if the Reynolds number is large, one must in each approximation solve a partial differential equation, which in the first approximation is Prandtl’s boundary-layer equation, but with unconventional boundary conditions. It is again convenient to work with the stream function rather than the velocity components and to eliminate the pressure from the Navier-Stokes equations by cross-differentiation. Then the flow is governed by the vorticity equation; in dimensionless form

+

(2.27)

Here the coordinates are again referred to a characteristic transverse dimension, and the stream function to the volumetric flow rate per unit thickness in half the channel. The Reynolds number R equals that flow rate divided by the kinematic viscosity. Then the boundary conditions are* *=+1,

+y

=0

at

y = EX).

(2.28)

As before, we contract the abscissa x, replacing it with X = E X , which transfers the perturbation parameter E from the boundary conditions to the

* Here we must recognize that linear terms in X do not contribute to the stress. * W e could require that #x vanish at the wall instead of # y , but that would complicate the calculations.

Milton Van Dyke

12

differential equation: QYYYY

+ 2E2Qxxyy+ E4Qxxxx = ER

6 = *l,

Qy = 0

at y

=

*f(X).

(2.29)

a. Reynolds Number of Order Unity If the Reynolds number remains fixed as E + 0, the differential equation (2.29) reduces to simply Qy;yy = 0, and the first approximation is then

lp)= i [ y / f ( X ) ]- i ( y 3 / f 3 >= ;7)

-

5q3.

(2.30)

This gives at any station the quasi-cylindrical approximation for the longitudinal velocity-the parabolic Poiseuille profile appropriate to flow with the prescribed flux rate through a parallel channel having the local width of that station. Again, however, because we have worked with the auxiliary function Q, which satisfies the continuity equation, this also yields the correct first approximation for the small transverse velocity component. To proceed systematically to higher approximations, we expand the stream function in powers of E . Equation (2.29) shows that, except in creeping flow, when the Reynolds number is negligibly small, odd as well as even powers of E now appear in the expansion (2.31)

Then substituting into the full problem (2.29) and equating like powers of E yields the sequence of problems

Blasius (1910) calculated the second approximation and used it to predict the separation in an exponentially growing channel. (He actually worked with the “primitive” variables u, 0,and p rather than Q and did not explicitly introduce the small parameter E, except in explaining the difference between what we have called slight and slow variations.) Abramowitz (1949) considered a special case of the third approximation. The channel-wall slope and the Reynolds number appear in the second approximation, and the

Slow Variations in Continuum Mechanics

13

wall curvature and the square of the Reynolds number in the third, giving

+ = -21( 3 7

-

3 v 3 )+ -~Rf'(577 280

- llv3

+ 77'

-

7')

+ 4488q7+ 1 1 5 5 ~ 987") ~ - T ( 1 2 1 3 ~- 3 2 7 9 + ~ 323477' ~ - 1518v7+ 3 8 5 ~ 357")] ~

x [f2(28757 - 8 2 2 2 ~ ~ 87787'

-

-

By delegating the mounting labor to a computer, Lucas (1972) calculated the 13th approximation for the general shape illustrated in Fig. 1 and additional terms for special shapes. He verified that the pattern apparent in (2.33) continues: the term + ( " + ' ) that multiplies E " contains alternate powers of R u p to R", odd powers of 77 up to 774"+3, and all possible n-tuple products o f f and its first n derivatives that involve just n differentiations. This last pattern is encountered in the potential solution (2.17), where only even powers of E appear. Lucas has calculated 25 terms of the series for a wedge with f ( X ) = X and tested it against the exact solution of Jeffery and Hamel (Rosenhead, 1963, p. 144). That solution in terms of elliptic functions has an infinite number of branches; but the Blasius series singles out the "principal" branch that includes plane Poiseuille flow and describes symmetric profiles with at most a single region of flow reversal near each wall. When the Reynolds number is negligibly small, the Blasius series (2.33) reduces to an expansion in powers of E ' . For the wedge, Lucas finds that it converges for E < 1. He extends its convergence by recasting it in powers of the semivertex angle tan-' E . By analyzing the coefficients, he finds that convergence is then limited by a square-root singularity at a semivertex angle of 129". This is the angle at which Fraenkel (1973) showed that the Jeffery-Hamel solution becomes singular, its derivative with respect to angle being infinite. When the Reynolds number is large compared with unity (say greater than 10 for practical purposes), the Blasius series reduces to an expansion in powers of ER. For the wedge, Lucas finds that it converges beyond the onset of separation (at E R = 4.712, according to Fraenkel), until it reaches a square-root singularity (at E R = 5.461). For other channels the convergence is not easy to judge. Among various shapes, Lucas shows the skin friction along a channel of sinusoidally varying width (Fig. 4) at high Reynolds number (that is, in the limit R -+ co, E + 0, with E R fixed). For ER = 4, the flow does not separate, and successive approximations appear to converge when recast into Pad6 approximants. For E R = 8, however, there appears to be a considerable extent of separated

14

Milton Van Dyke

flow; and Fig. 4 shows that the Pade approximants formed from 7, 9, and 11 terms of the series have not yet settled down (if indeed they ever will). b. High Reynolds Number Suppose that the Reynolds number is so large that the product E R is not small. We have seen that the Blasius series may converge; if so, an accurate solution can be found from sufficiently many terms. To attack the nonlinear problem directly, however, we must retain terms of order E R in the problem (2.29). With an error of order E * it becomes (1)

= ER($y

L ; : $ $(I)

=

*I,

(1)

$x,y

(1)

- $x

$y=0

(1) $yyy)r

at y

=

*tf(X).

(2.34)

This approximation was first proposed by Williams (1963), who observed that the differential equation is that of Prandtl's boundary-layer theory. However, it satisfies zero-velocity boundary conditions at both walls rather than being matched to an outer inviscid flow, and the pressure (which was here eliminated by cross-differentiation) is not imposed, but is found in the course of solution. Williams searched for self-similar solutions in symmetric channels and found only the slow-variation version of flow in a wedge. Later Blottner (1977) used a marching finite-difference scheme to calculate non-self-similar solutions of (2.34) for various channels. In particular, he treated a family of symmetric channels expanding slowly from a parallel upstream section through a cosine variation to a parallel downstream section of double the upstream width. For the shortest of his transitions he found separation at the walls followed by a short region of reversed flow. Similarly, Eagles and

Cf

0

FIG. 4. Skin friction in periodic channel y approximants formed from 2 N + 1 terms of series.

=

* 2 sin E X for E R = 8. [ N / N ] Pad6

Slow Variations in Continuum Mechanics

15

Smith (1980), as a basis for studying the stability of nonparallel flows, solved numerically the flow through a channel expanding to three times its upstream width according to + y = 1 + tanh( E X ) . Their forward-marching scheme showed separation starting at ex = 0.35 when E R reaches 15.5, and extending over 0.2 < E X < 1.5 when e R = 18. B. AXISYMMETRIC SHAPES Axisymmetric problems of slow variation can be treated in exactly the same way as symmetric plane problems, and the results have a remarkably similar structure. We illustrate this by examining the solutions of three different problems. Our region of interest (Fig. 5) is the axisymmetric counterpart of the strip considered previously. We describe it by using the same slowly varying function f(e x ) , but replace the transverse Cartesian coordinate y by the cylindrical radius r. 1. Potential Flow through Pipe of Varying Radius

Of the various physical interpretations that could be given to the problem of the potential in a plane strip, that of temperature distribution has no realistic counterpart here, and so we consider potential flow. Then it might seem natural to work with the velocity potential, but we have seen that we gain half a step by using the stream function $ instead. Then the velocity components are given in cylindrical coordinates by u = (l/r)(d$/dr),

u = -(l/r)(d$/dx).

(2.35)

In axisymmetric flow the stream function no longer satisfies the Laplace equation, but instead the equation D2$ =

ILrr

- (W r ) +

In terms of the contracted abscissa X

=

$xx

= 0.

(2.36)

x, this becomes

r ( $ r / r ) r = -&xx.

FIG.5 . Slowly varying axisyrnrnetric shape.

(2.37)

16

Milton Van Dyke

Here we have combined the two r-derivatives in (2.36) to make it evident that again only quadratures are required at any stage of the approximation. As boundary conditions we set = 0 on the axis r = 0 and, in appropriate dimensionless variables, CC, = on the wall of the pipe r = f ( X ) . Then the first approximation is seen to be simply @ = r2/2f2(X), which gives the axial and transverse velocities

+

u = 1/f2(X),

2,

(2.38)

= E(Tf/f3).

Again the first of these is the hydraulic approximation of constant axial velocity at each cross section, and the second is our bonus obtained by working with Carrying the approximation to third order gives

+.

+ = (r2/2f2) +

&E’(

l/f2)”r2(f2 r2)

+ ~ ~ { ~ [ f ’ ( l / ” ) ’ ’ l ’ ’-r r~2() f-~& ( ~ / f ~ ) ( ~ b-~r4)} (f+ ~ o(E~). (2.39)

This expansion has a structure similar to that of (2.16) for the corresponding plane flow. The resemblance is increased by recasting into Lucas’s form, carrying out the indicated differentiations, and replacing the radius r by the fractional distance 7) = r/f to the wall. This gives $=L

(3y2- f f ” ) ( v 2 v4) + E 4 x [&( 18y4- 32ff”f” + 7f2T2 + 8f2fY- f 3 f ‘ 4 ’ ) (

2 7

-

2 + ’ 2

sE

&(,Of4

- 72ff”f”

+ 9f2Y2+ 12f2fY

-

v 2 - v4)

f3f‘”)(77’

-

v6)]. (2.40)

Here + / q has exactly the same form as I) itself in the corresponding solution (2.17) for plane flow, with only the coefficients being different. We can check this result and its convergence in the case of a conical pipe. The exact solution is given by a point source at the vertex, so that for a cone of semivertex angle a = tan-’ E it is 1 1-COS~ 1 1,21 - ( 1 + E2r2/X2)p’/2 *== $1 + E 1 ( 1 + &2)1/2 - 1 . (2.41) 2 1 - cos ff Expanding this for small E reproduces the result of setting f ( X ) = X in (2.40). This shows that the series converges for E < 1 in this case. 2. Torsion of Shaft of Varying Radius Let the axisymmetric shape of Fig. 5 now represent a solid elastic shaft that is twisted by couples applied at its distant ends. According to the principle of Saint-Venant, the precise way in which they are applied will significantly affect the deformation only within “boundary layers” that extend a couple of radii from the ends. Elsewhere each particle moves only

Slow Variations in Continuum Mechanics

17

in the tangential direction during twist (Timoshenko and Goodier, 1951, Section 104). Then the two nonzero components of stress can be expressed in terms of Michell's stress function 4 according to Tre =

-4

x 1 r2,

=

(2.42)

+r/r2.

The stress function satisfies the differential equation (2.43) 4rr- 3 ( 4 J r ) + 4xx= 0. The boundary conditions require that 4 increase by a constant amount

from the axis to the surface, say in dimensionless form at at

4=0

4=1

r=0, r=f(Ex).

(2.44)

Proceeding just as before yields for the third approximation

4

+

= ( r 4 / f 4 ) &E'(

1 / f 4 ) " r 4 ( f 2- r 2 )

+ ~ " ( ~ [ f ( l / f " ) " ] " r " -( f r' 2 ) - L(1/f")'"'r"(f4 - r")}. 384

(2.45)

In terms of the fractional radius 77, this becomes

4

=

+ i E 2 ( y-jy)7744(1 2 + E4[~(i00y4 - 112f2y+ 1 3 f 2 y 2 + 16f2f'f"' - f 3 ~ ( ~ 9 ~ 4 ( 1 - 772)

774

772)

-

&(210y4- 180f2y

+ O(E6).

+ 15f2y2 + 20f2f'f"'

+

- f 3 f ( 4 ) ) 7 7 4 ( 1 - q4)1 (2.46)

Here 4/773has exactly the same form as the of (2.17) for plane potential flow or the + / T of (2.40) for axisymmetric potential flow, but with different coefficients. This similarity is not surprising, because all three problems are members of the family of "generalized axisymmetric potential theory," of the form at r = 0, (2.47) at r = EX), with N being 0, 1, and 3, respectively. The approximate solution of this generic problem for slow variations is, to second order,

4

=

N+l [ ( N + 2)f" - f f " ] 7 7 N + 1 ( 1 - T ~ ) . (2.48) 2( N 3 )

$J+'+ & 2

+

Foppl found the exact solution for torsion of a conical shaft (Timoshenko and Goodier, 1951, Section 104) as

'

2 - 3 cos 6 + C O S ~0 = 2 - 3 cos a + cos3 a 2 - 3(1 E2r2/X2)-1/2 (1 E ~ ~ ~ / X ~ ) - ~ / ~ (2.49) 2 - 3(1 + E ~ ) - " + ~ (1 + E ~ ) - ~ / ~

+

+ +

18

Milton Van Dyke

Expanding for small E checks (2.46) whenf(X) the series converges for E < 1 .

=

X and shows that again

3 . Laminar Flow in Tube Blasius concluded his pioneering paper of 1910, which is devoted mainly to plane flow, by treating the corresponding axisymmetric problem of laminar flow through a straight pipe of slowly varying circular cross section. It is again convenient to work with the vorticity equation. In cylindrical polar coordinates and dimensionless variables, the Stokes stream function $ satisfies

(2.50) where the elliptic operator D’ is given by (2.36). The boundary conditions are =

as at

O(r’)

*=’’,

I+!I~

=O

r + 0, r=f(~x).

(2.51)

Again contracting the axial coordinate to X = E X transfers the perturbation parameter E from the boundary conditions to the differential equation, which becomes

x

(a’*

ar2

a*

r ar

+ E’%).

ax

(2.52)

a. Reynolds Number of Order Unity If the Reynolds number remains fixed as F + 0, the differential equation (2.52) reduces to (a*/ar’ - r-’a/ar)’+ = 0. This gives again at each station the local Poiseuille flow described by

9 = [ r ’ / f ’ ( ~ ) -] ;[r4/f4(x)] = 7’

- $17“.

(2.53)

Blasius (1910) found the second approximation, and Manton (1971) has calculated the third. Kaimal (1979) has generalized Manton’s results to a dilute suspension of solid particles and in so doing has tacitly corrected an error in one of Manton’s coefficients.* (Unfortunately, he has at the same time introduced typographical errors into two of the other coefficients.) We *The expression below Manton’s Eq. (4.7) should be a , = ‘f&

+ 13Q).

Slow Variations in Continuum Mechanics

19

may write Manton’s corrected result as $=

1

(

’

v2- -v4 + -~ERf-(477’ :

)

+ ~’{:(5f’

1

8

f

-9v4

+

- - ( 1 5 6 ~~ 4 3 5 ~ 450$ ~

f

””

-ff”)(~’- 2q4 + q6)+ - -(81877’ R 2 f’ 5400

+ 2 6 5 0 -~ 1~ 4 2 5 +~ 39017”~ f”

+ 6776 - 7’)

387”) -

22577’

+ 6077”

- 2395q4

- 677”)]}

+ O(E~). (2.54)

Chow and Soda (1972) independently undertook exactly the same analysis, but their result for the third-order terms that are multiplied by R bears little resemblance to that of Manton (which we have checked and corrected); in particular, they do not include the effect of the longitudinal curvature, as represented by f”. This form of Manton’s expansion is the axisymmetric counterpart of Lucas’s solution in the form (2.33). Comparison shows that the pattern is remarkably similar in the two cases, especially for the terms independent of the Reynolds number. Aside from an extra factor of 77 in the axisymmetric case, the plane and axisymmetric expansions for creeping flow at R = 0 have identical structure with respect to the perturbation parameter E, the radius f and its derivatives, and the powers of the fractional radius 7.The terms involving the Reynolds number have a different structure, however. In fact, the axisymmetric form is obtained from that for plane flow by replacing R throughout by R / f and then the two expansions differ only in their coefficients. b. High Reynolds Number When the Reynolds number is so great that the product ER is of order unity, the largest nonlinear terms on the right-hand side of (2.52) must be retained in the first approximation. Again this gives the equation of classical boundary-layer theory. Williams (1963) studied that approximation by replacing the radius with the fractional radius 77 = r / f ( EX) and discovered a reduction to an ordinary differential equation when the tube radius f ( EX) is an exponential. Thus self-similar solutions exist for laminar flow through a tube whose radius varies slowly and exponentially. Williams solved the resulting ordinary differential equation numerically and plotted the family of velocity profiles for exponentially converging and diverging tubes, showing reverse flow near the wall in the more divergent cases. Daniels and Eagles (1979) investigated the problem in greater detail and discovered multiple solutions. However, they believe that only the solutions on the branch computed by Williams, which includes Poiseuille

20

Milton Van Dyke

flow in a straight tube, are physically realistic. That is the branch represented by Manton's series. This self-similarity can be regarded as the axisymmetric counterpart of the Jeffery-Hamel similitude for plane flow, except that it is only a first approximation for a slowly varying exponential pipe, whereas the JefferyHamel solution is exact for any wedge angle. Nor has the solution been found in closed form, but only by numerical integration. This limited similitude can be seen in Manton's series (2.54). When f ( X ) = exp(aX) the factor ERf '/f in the second term becomes independent of X and equal to Daniels and Eagle's similarity parameter y = EAR. Likewise the factor & ' R 2 Y / f , as well as E 2 R 2 f t 2 / f 2in, the third term becomes the square of that parameter; and higher powers would appear in subsequent terms. On the other hand, the terms in E~ that do not contain R represent the start of a correction to the self-similar solution, of order ~~u~exp(2aX), which Daniels and Eagles discuss briefly. For a tube with other than exponential variation, nonlinear partial differential equations would have to be solved numerically. This has not yet been done, as it has been for plane flow. Instead, Eagles and Muwezwa (1986) have proposed expanding in powers of A = E R as well as E'. This evidently yields just a reordering of the expansion (2.54) of Blasius and Manton. Eagles and Muwezwa actually treat only the leading term of the expansion in powers of E ~ so , that they are solving the classical boundarylayer equations by expanding in powers of ER. Thus their first two terms are those in (2.54), and their third term is the part of (2.54) proportional to s 2 R 2(which confirms our correction of Manton's solution). They break new ground by giving also the term in E ~ RExamining ~ . the coefficients of this four-term series for a divergent and a convergent-divergent tube, they suggest that the series has a finite radius of convergence rather than being merely asymptotic, and they use it to compute velocity profiles up to E R = 10. C . ANTISYMMETRIC PLANESTRIP

So far we have considered only shapes having a straight central axis of symmetry. We turn now to unsymmetric shapes, beginning with the simplest case of a plane strip having a curved centerline. In general, both the width and the direction of the strip will vary slowly along its length (Fig. 6a). However, having already treated varying width for a straight centerline, we can illustrate the new features associated with curvature of the centerline by restricting attention to a strip of constant width (Fig. 6b). We could then easily treat the general case by combining the techniques used for the two special cases. Cartesian coordinates are no longer appropriate to the problem, for although the strip changes direction only slowly, it may eventually reverse

Slow Variations in Continuum Mechanics

21

FIG. 6. (a) General case of slowly varying strip with curved centerline; (b) special case of constant width.

its course. Instead, we employ as coordinates the curvilinear distance s along the centerline and the distance n normal to it (Fig. 6b). Then we describe the shape of the centerline implicitly by giving its curvature K (or radius of curvature 1 / ~ as ) a function of s. We choose to reckon the curvature positive if the center of curvature lies at negative n [which is the convention used in boundary-layer theory (Rosenhead, 1963, p. 202), but the opposite of that used by Wang (1980a) and Van Dyke (1983)], and we choose the length scale such that the width of the strip is 2.

1. Potential in a Meandering Strip Again we seek first a solution rC, of the Laplace equation. We require it to assume the values k1 on the upper and lower boundaries. The line element is ( 1 + Kn)ds2 + dn2, so that Laplace’s equation in our curvilinear coordinates system becomes v2*

=-

(2.55)

If we assume merely that the curvature K is small, we are again making the approximation of slight variations. The result does not fail for large s, as did the corresponding approximation for the symmetric strip, because for our special case of a strip of constant width the coordinates conform to the boundaries and we do not have to transfer the boundary conditions. However, this approach still has the disadvantage that beyond the first approximation we must at each stage solve a Poisson equation, which can be done only for a specific shape. This difficulty disappears if we assume that the curvature is not only small, but also slowly varying. We express this by introducing the contracted

22

Milton Van Dyke

longitudinal coordinate S = E S and setting K(S)

=

Eg(s).

(2.56)

Then the Laplace equation (2.55) gives

a

-(1+ an

a+

Egn)-+ an

E

a

1

a+

as1 + EgnaS

= 0.

(2.57)

This evidently reduces to (Clnn = 0 in the first approximation and to the same with a known right-hand side at each subsequent step. Thus again only quadratures are required to treat a general shape. In contrast to the previous problem for a symmetric strip, the expansion here evidently proceeds in odd as well as even powers of E. We readily calculate the fifth approximation as

+ = n + $&g(l- n2) - fs2g2(n - n’) - &e3[(4g3 - 6g”)(l - n2) + ( g ” - 6g3)(1 - n‘)] + &‘[&(4g4 - 7g” - 13gg”)(n - n’) + &7gf2 + llgg” - 24g4)(n - n’)] + O ( E ’ ) .

(2.58)

Here the N t h term, proportional to E ~ - ’ , contains alternate powers of n up to the Nth; and (except in the first term) each of these is multiplied by a linear combination of all products of the function g and its derivatives that contain N - 1 - 2i factors and 2i derivatives, where i = 0, 1 , . . .. If the function g ( S ) is taken to be a constant, the strip becomes a thin annulus, with g = 1 in the notation of Fig. 7. Then the exact solution is

FIG. 7. Thin annulus as a meandering strip.

Slow Variations in Continuum Mechanics

23

given by a vortex at the center as $=1-2

1~n)] In[(l + ~ ) / ( + In[(l + & ) / ( I - E ) ] '

The expansion of this result in powers of case and converges for E < 1.

E

(2.59)

checks our series (2.58) in this

2. Laminar Flow in a Meandering Channel Consider now the steady plane laminar flow through a channel having the shape of Fig. 6b. The Navier-Stokes equations yield for the stream function $ the vorticity equation (2.60) where the Laplacian V2 is given by (2.55). Here the Reynolds number R is again the volumetric flow rate per unit thickness in half the channel divided by the kinematic viscosity and the stream function is normalized to $ = *l at the channel walls n = *l. Wang (1980a) has treated this problem by assuming that the curvature is small, with K = Ek(s). In order to approximate for small E, he finds it necessary to assume also that the Reynolds number R is small, of order E. Then at each stage beyond the first he must solve a nonhomogeneous biharmonic equation. For the periodic channel with centerline curvature K = E cos As, he is able to proceed as far as the nonperiodic terms in the third approximation. a. Reynolds Number of Order Unity

As with the potential flow of the previous section, the process is greatly simplified by assuming that the curvature is not only small, but also slowly varying, as in (2.56). If the Reynolds number remains fixed as E + 0, the vorticity equation (2.60) reduces to simply a4$/dn4 = 0 at the first stage and to its nonhomogeneous counterpart in higher approximations, so that again only quadratures of polynomials are required. The fourth approximation is readily found (Van Dyke, 1983) as 1 1 $ = -(3n - n') + - E g ( S ) ( 1 2 4

1

x (67 - 2n2 - n") (1 - n2)2- E' 1 + -(189n 3150

- 13n' - 5n5)Rgg'

x (4663 - 2940n2

R

- n2)' - E'

1 + 240 + 3n2)g"

- 33n2)g3 -(1

1 + 1,478,400

1

+ 842n4 + 4n6 - 9ns)R2g" (1 - n2)' + O ( E ~ )(2.61) .

24

Milton Van Dyke

For creeping flow, at negligible Reynolds number, this series has the same structure as (2.58) for potential flow, except that each polynomial in n is of two-higher degree in order to accommodate the no-slip condition at the wall. For an annulus the curvature K = Eg is constant, and the exact solution has the form $ = A(1+ Egn)2 ln(1

+ Egn) + B ( 1 + Egn)2 + C ln(1 + Egn) + D.

(2.62)

Here the coefficients A, B, C, D are rather complicated functions of Eg that are determined b y imposing the boundary conditions (2.28) that $ = f 1 and $n = 0 at n = *l. Expanding for small E checks our series (2.61) and shows that, as in the potential problem, it converges for E < 1 in this case. b. High Reynolds Number When the Reynolds number is so great that the product E R is of order unity, the leading nonlinear terms in (2.60) must be retained in the first approximation. With an error of O( E ) it becomes $nnnn

=

ER($n$snn

-

$s$nnn)*

(2.63)

This is again the classical boundary-layer equation. It is identical with the equation (2.34) for a symmetrical channel (with X replaced by S and y by n ) , because longitudinal curvature has no first-order effect in Prandtl’s boundary-layer theory. Therefore the self-similar solutions of Williams (1963) for narrow wedges and the non-self-similar solutions of Blottner (1977) and Eagles and Smith (1980) for symmetrical channels can be applied to channels with slowly curving centerlines. The only change is that neglect of centrifugal pressure then imposes an error of O ( E ) compared with O( E ’ ) in the symmetric channel. 111. Three-Dimensional Slender Shapes

The extra freedom that arises from allowing extension into a third dimension enormously increases the variety of possible slowly varying shapes. Only a few special configurations have yet been calculated, and those mostly for laminar flow. Much further study of three-dimensional slowly varying shapes is to be expected. We first direct attention to slender shapes, those whose two transverse dimensions are both small compared with their longitudinal dimension. The centerline can meander through three dimensions so that (as for a helix) it has torsion or twist as well as curvature, and the cross section can vary slowly in shape, proportion, size, and orientation. In unpublished work

Slow Variations in Continuum Mechanics

25

Todd (1978, 1979, 1980) has studied steady laminar flow through a slender pipe of such general slowly varying shape. He derives the governing equations and then describes the scheme of successive approximations that yield in turn the basic Poiseuille flow in the axial direction, the first approximation to the secondary motion in the cross section, and the first correction to the axial flow. He devotes special attention to the axial pressure gradient and is able to draw various conclusions without further specifying the shape of the pipe. To obtain more definite results, however, he considers special cases, mostly with circular or elliptic cross sections.

A.

STRAIGHT CENTERLINE

We choose to classify slender shapes first according to whether the centerline is straight. This is an essential distinction because, as we have already seen for plane strips, we can often use conventional coordinates if the centerline is straight, but most otherwise use curvilinear (and in general nonorthogonal) coordinates. 1 . Flow through Pipe of Varying Elliptic Section

Our first example of a nonplanar solution was the potential flow (2.38) through a straight pipe of slowly varying circular cross section. We now consider the generalization to a pipe of elliptic cross section, whose major and minor axes may vary slowly and independently. Thus it is described, in terms of the contracted axial coordinate X = E X , by y ’ / a 2 ( X )+ z’/b’(X)

=

1.

(3.1)

a. Potential Flow We used the Stokes stream function for the circular pipe, but no single stream function exists for a three-dimensional flow. We therefore treat potential flow through the elliptic pipe by using the velocity potential 4. It satisfies Laplace’s equation in the form +y,,

+ 4 z z + E 2 4 x x = 0,

(3.2)

and the flux through the pipe is prescribed, in dimensionless variables, as

The hydraulic approximation, with uniform axial velocity across each section, gives u =

c $ ~=

l/a(X)b(X).

(3.4)

26

Milton Van Dyke

The second approximation for the velocity potential was calculated by Olson (1971), and his analysis is reported by Sobey (1976), who uses it as the basis for treating inviscid flow with slight shear. With an error in sign corrected in Sobey’s expression for it gives the three Cartesian velocity components 1 1 ab 8 a’ E2Y a b

u = - - - E’ =

[(-$)

’ ( a 2- 4y2) +

b‘ + o ( & ~ ) , = E2Zab +

(3.5) 0(&3).

Rewritten in terms of the fractional transverse coordinates 7 = y / a ( X ) and 4‘ = z / b ( X ) ,these are

=

[

a2(L)’(1 ab 8 a2b a’ &-? + 0(&3), ab

= -- - & 2

- 477’) =

+ b2($)’(1

+O(E~),

- 412)]

(3.6)

b’ ab

&-l+ o ( & ~ ) .

It is not clear that the next approximation will have the simple form of the first and second, a polynomial in y and z with coefficients depending on X . We see that having to work with the velocity potential has left us half a step behind working with the stream function, where the second approximation provides the secondary terms in the transverse as well as the axial velocities. However, we are ahead of working with the primitive variables, as we must in the next example. b. Laminar Flow The corresponding laminar flow has been calculated by Wild ef a/.(1977). Since neither a velocity potential nor a single stream function exists, they work with the Navier-Stokes equations for the three Cartesian velocity components and the pressure. The first approximation for the axial velocity u is Poiseuille flow, which has a simple closed-form solution for an elliptic pipe. In terms of the fractional coordinates 77 and 5, it is u = ( 2 / a b ) ( l- 9’-

6’).

(3.7)

Next, using the continuity equation and eliminating the pressure between the y - and z-momentum equations by cross-differentiation yields the first approximations for the transverse velocity components: 2a‘ u = &--77(1 -7 ab

2 -

52),

2b’ ab

w = s-6(1

-

77

2

-

62 ).

(3.8)

These exhibit a kinship with their counterparts (3.6) in potential flow.

Slow Variations in Continuum Mechanics

27

In the third step, a correction to u of order E R is calculated from the axial momentum equation. This gives L

u = - ( l - v 2 - 5 2) ab x [ I + E R ( C ~+ c , q 2 + c2v4+ c,c2 + c414 + c5v2l2) + O ( E ~E ,~ R ~ ) ] . (3.9)

We refer the reader to Wild, Pedley, and Riley’s equations (2.13) and (2.14) for the coefficients co through c 5 , which depend on the contracted axial coordinate X in fairly complicated fashion through the functions a ( X ) and b(W. 2. Flow through a Twisted Elliptic Pipe

Todd (1977) has treated the fully developed laminar flow through a straight pipe whose cross section is constant, but rotates slowly about its centroid in the axial direction. He considers a general cross section, but derives detailed results only for the ellipse. That cross section is described by y 2 / a 2 z 2 / b2 = 1, where x, y, z form a slightly nonorthogonal coordinate system that twists with the pipe. As shown in Fig. 8, those coordinates are

+

FIG. 8. Coordinate systems for slightly twisted elliptic pipe.

Milton Van Dyke

28

related to the Cartesian system 3, j j , Z by x = 2,

y

=y

cos 8 + Z sin 8,

z

= Z cos 8

-

j sin 8.

(3.10)

Here the angle 8, shown in Fig. 8 for right-handed twist, is a slowly varying function of the axial distance x. We make that slow variation explicit by setting

x=

8(x) = @ ( X ) ,

(3.1 1 )

EX.

Here 0 is a normally varying function of its argument and is simply X for a constant twist angle E. a. Potential Flow Again we examine first the simpler problem of potential flow. The velocity potential 4 satisfies the Laplace equation, which in our twisting transverse coordinates y, z and contracted axial coordinate X = E X becomes 4yy

+ 4 z z + E 2 [ 4 X X + 2@’(z4xy- Y 4 z x ) + @”(z4, - Y 4 Z ) + @’2(z24yy - 2YZ4yz + Y24*z- Y 4 y - z4,)l = 0.

(3.12)

The condition of tangent flow at the wall of the pipe is b2y4,

+ u ’ z ~ ,= E’@’(u’ - b2)yz[4x + @’(z+, at

(y’/la2)

-

y4,)]

+ (z’/b2) = 1 ,

(3.13)

and the condition of prescribed flux through the pipe is, in normalized form, (3.14) If the twist is negligible, the solution is simply a uniform parallel flow with 4 = x = X / E .Then the governing equations [(3.12)-(3.14)] show that for small twist angles the perturbation expansion proceeds in powers of E ~ Carrying out the second and third approximations gives

4

=

U’ - b2 E @ ’ ( X ) ,yz P + &3[(4@’3- (3’”) a2+ b x (c,a2yz + c2y3z+ c3yz3)+ Or@” x (c4a4 c5a2y2 + C6U2Z2 + c,y4 + c*y’z’ + c9z“)] + * * .

X -+ E

+

-

(3.15)

Here c,, c2, . . . , c, are algebraic functions of the thickness ratio p = b / a of the cross section. The first two terms of this expansion show that the streamlines twist more slowly than the pipe itself, in the ratio ( a ’ - b2)/(a2+ b2) = ( 1 - p 2 ) / ( 1 + p’). This factor agrees with one’s physical perception that the streamlines will not twist at all through a round pipe ( p = l ) , but must twist nearly as fast as the elliptic section when it becomes flat ( p << 1 ) .

.

Slow Variations in Continuum Mechanics

29

b. Laminar Flow Todd has treated the corresponding problem for laminar flow by using the Navier-Stokes equations for the three velocity components and the pressure. He shows that those equations assume rather complicated forms in the twisted coordinate system x, y, z of Fig. 8. The basic solution, for Poiseuille flow through an untwisted elliptic pipe, gives with our normalization (3.14) the axial velocity ~ = 2 [ -y’/a’-~’/b’]. 1

(3.16)

Twist of the pipe induces a crossflow. Todd finds that for small twist the velocity components in the directions of the major and minor axes of the cross section are given by 2,

=

EO‘(X)

3

+ 2p2 + 3p4

( 1 - P2)(3+ P 2 ) w = &O‘(X) 3 2p2 3p4 yu,

(3.17)

+

+

where u is given by (3.16). These results show that the streamlines again twist more slowly than the pipe itself, now in the ratio (1 -

3

+ 2p + 3p2 +

2p2

+ 3p4‘

(3.18)

This is smaller than the ratio for potential flow; evidently the no-slip condition at the wall inhibits the spiraling of streamlines. Powers of the Reynolds number enter into high approximations, which contain (in contrast to the potential flow) both even and odd powers of the twist parameter E. Thus Todd shows that the correction to the basic axial velocity (3.16) is in general proportional to FRO’;however, in the special case of the ellipse that secondary term vanishes, leaving only tertiary correction terms proportional to E ~ O E’O’’, ‘ ~ , E ~ R ’ O ’ ’and , E’R‘O’’ to be added to (3.16). These are multiplied by polynomials in y and z of degree as high as 10. The crossflow velocities of (3.17) have secondary corrections proportional to E’RO’’ and E’RO”, times polynomials in y and z. The coefficients of all these polynomials are rational functions of the thickness ratio 0, which Todd finds much too complicated to calculate except numerically by computer.

3. How through a Spiraling Circular Pipe Kotorynski (1979) discusses the steady laminar flow through any slowly varying three-dimensional channel and then treats as a specific example the slowly spiraling circular pipe shown in Fig. 9. This pipe can be viewed as being generated by translating a circle of unit radius so that its center moves along a helix of radius a while its plane remains normal to the axis of the

30

Milton Van Dyke

FIG.9. (a) Slowly spiraling pipe and (b) associated coordinates. Shown dotted is the nearly circular straight pipe that is tangent at station xo.

helix. The cross section of the pipe is then not quite circular in a plane normal to its own helical axis, which is what is usually understood as a coiled circular pipe, but the difference affects only third-order terms. As in the preceding example, we introduce a slightly nonorthogonal coordinate system x, y, z that spirals with the pipe, as shown in Fig. 9b. It is related to the Cartesian system 2, j , i by x = 2,

+

y = jj cos( E Z ) ,? sin( E X ) - a, z = 2 cos( &a)- j j sin( & a ) .

(3.19)

a. Potential Flow Again we consider first the simpler problem of potential flow. In the spiraling coordinate system just introduced, with the axial coordinate contracted to X = EX, Laplace's equation for the velocity potential becomes 4yy

+ 422 + E"4XX + 3 2 4 , - (Y + aMxJ + z24yy- 2(Y + a)z&,, (3.20) + (Y + a)*4,, - ( Y + a M y - 2423 = 0.

The condition of tangent flow at the wall of the pipe is y&

+ z4, = a z [ ~ +$ ~z4,, - ( y + a ) h ]

at y 2 + z 2 = 1 , (3.21)

and we prescribe the flux through the pipe in dimensionless form as

11

4,dydz

=E

11

4 x d y d z = T.

(3.22)

Slow Variations in Continuum Mechanics

31

When E is negligibly small, the solution is a uniform stream along the axis of the helix, with 4 = x = X / E . Higher approximations proceed in powers of E’, and one easily calculates the third approximation as 1

=-X E

+ E U Z + -81E ~ U [ Y ~ +Z z 3

-

~ U Y Z-

( 3 + 8 a 2 ) 2 ]+ . . .. (3.23)

Here the first two terms give the velocity components to order E as u = 1 , u = 0, w = &a.This represents at each axial station a stream that is uniform and parallel, with unit speed, but is slightly inclined to the x axis, at an angle w / u = &a. We see from Fig. 9 that this is just a uniform flow directed along the spiraling centerline of the pipe, rather than along its straight central axis. b. Laminar Flow Kotorynski treats the corresponding laminar flow by using coordinates that differ slightly from ours in orientation. (His y and z have the same spiraling origin as ours, but remain parallel to themselves rather than rotating.) His basic solution is, of course, the Poiseuille flow through a straight circular pipe; with our notation and normalization this gives u=2(1-r2)=2(1-y2-z2).

(3.24)

After some analysis Kotorynski finds the first approximation to the transverse elocities; in our notation again (and with a missing factor a restored) these are u=o,

w=2m(l-y

2

2

- 2 ) .

(3.25)

We notice that this flow field again represents at each axial station simply a parallel stream slightly inclined to the central axis. Rather than a uniform stream, it is now a Poiseuille flow with parabolic profile, but it is again inclined at an angle w / u = E U . Hence the solution to this order is nothing more than the basic Poiseuille flow in a straight circular pipe, but one that is tangent at each axial station to the actual curved pipe, as indicated in Fig. 9a. True effects of curvature would appear in subsequent approximations (together with the effects of the slight deviation from circularity of the tangent pipe). It seems likely that the next approximation could be calculated without great labor. Even in this problem, however, where the axis of the pipe never departs far from the central x axis, we would probably simplify the analysis b y recognizing that the centerline is curved and using the curvilinear distance s along it as the axial coordinate (as we did for plane flows in Section 11, C and will do for three-dimensional shapes in the following sections).

32

Milton Van Dyke B. CURVEDCENTERLINE

We turn now to slender three-dimensional shapes with curved centerlines. As in the corresponding plane problem, we restrict attention to unchanging cross sections, and it appears that in fact only such shapes have so far been treated in the literature. 1. Twist of an Elastic Ring Sector Gohner (1930, 1931a) has considered the elastic properties of a helical spring in tension or compression. He shows that if the pitch of the helix is small, it can be neglected with only very slight error. Thus the spring is modeled by a sector of an annular ring whose ends are subjected to equal and opposite forces directed along the axis of rotational symmetry,.as shown for a circular cross section in Fig. 10. The stress field is the same in each cross section around the ring. If the spring is loosely coiled, the stresses are, to a first approximation, those due to torsion in a straight shaft of the same cross section. Gohner finds the effects of curvature by expanding in powers of the ratio E = a / L of a typical cross-sectional dimension a to the coiling radius L. He computes the corrections of order E and E’ for circular, elliptic, and rectangular cross sections. Timoshenko and Goodier (1951) present a concise and accessible account of the analysis for a circular section. The only nonzero components of stress are the shearing stresses T~~ and reyacting over the cross section. The equation of equilibrium in which they appear can be satisfied by introducing a stress function 4 such that G L ~d 4 ( L + x ) a~y ’

(3.26)

rex = ~-

where G is the modulus of rigidity. Substituting these into the two relevant

t

X

I F FIG. 10. Loosely coiled ring sector of circular cross section.

Slow Variations in Continuum Mechanics

33

compatibility equations and eliminating T~~ and T~~ by cross-differentiation yields the equation to be solved for the stress functions; in dimensionless form it is

$4

a2d, -+----ax2

dy

3.5

1+ Exax

2c

=

(3.27)

0.

Here c i s a constant that is at the end of the calculation evaluated in terms of the applied torque FL by integrating the shear stress over the cross section. Here our x differs in sign from Gohner’s (and from Timoshenko and Goodier’s equivalent 5) because we measure it away from, rather than toward, the axis of rotational symmetry. The condition that the shear stress be tangent to the boundary of the cross section is satisfied by requiring 4 to vanish on the boundary. When the coiling ratio E is negligible, the third term in (3.27) disappears, and the solution is 4 = ;c(l - x2 - y’), which is the stress function for torsion of a straight shaft. The effects of coiling are introduced by expanding in powers of E ; and thus the third approximation is found to be

4

= ;c(l - X’ - y2)[1

- $EX

+

&E’(X~

+ 5 y 2 - 15) + . .

a].

(3.28)

The shearing stress is greatest at the innermost point of the cross section, where it exceeds the value for a straight shaft by the factor 1 + dE + ;&2 - . . .

(3.29)

It would be interesting to calculate further approximations if only to check a remarkable closed-form expression for the greatest shearing stress, quoted by Timoshenko and Goodier (1951), that was “communicated to S. Timoshenko in a letter from 0. Gohner.” In a subsequent paper Gohner (1931b) treats a ring sector subject to a general loading at its ends. Timoshenko and Goodier summarize the analysis for the case of a circular cross section and for pure bending in the plane of the ring. 2. Laminar Flow through a Coiled Pipe Much attention has been devoted to flow through a coiled pipe because of its practical importance. The simple potential flow has apparently never been analyzed, but steady, fully developed laminar flow through a loosely coiled pipe has been calculated for several cross-sectional shapes as a perturbation of the Poiseuille flow through a straight pipe. Just as Gohner neglected the pitch of a helical spring, so most investigators have approximated the coiled pipe by a torus. In creeping motion, when the Reynolds number is negligibly small, there is no secondary motion, the flow being directed normal to the cross section. For a circular section Larrain and Bonilla (1970) expand the solution in

34

Milton Van Dyke

powers of the coiling ratio u = (1

- r2)[1 - ~

E

=

a / L to find the velocity as

cos 8 - &’(3

E T

- l l r ’ - lor2 cos 28)

+ - -1. 1

(3.30)

A quantity of interest is the flux ratio, the flow rate through the coiled pipe divided by that in a straight pipe under the same pressure gradient. This is given by Q c / Q s= 1 + Z E

+ . . ..

1 2 - 1 3 10248~

(3.31)

Thus the flux is actually increased by slight coiling because, as Larrain and Bonilla pointed out, the velocity distribution (3.30) shifts toward the inner wall of the pipe and so shortens the fluid path. When the Reynolds number is not small, Dean (1928) showed that to a good approximation the flow depends not upon the coiling ratio and Reynolds number separately, but only upon the combination K

= RE’/’ = ( 2 f i a / ~ ) ( a / L ) ’ / ~ .

(3.32)

This (or one or another related parameter) has become known as the Dean number. Here R = 2ua/ v is the Reynolds number based on the mean axial velocity U and pipe diameter. Dean worked instead with a fictitious Reynolds number based on the maximum velocity that would arise in a straight pipe under the same pressure gradient; when recast in terms of the actual Reynolds number, his approximation for the flux ratio is

Q c / Q s= 1 - 0 . 0 3 0 5 8 ( ~ ’ / 2 8 8 )+~ 0 . 0 0 8 1 9 ( ~ * / 2 8 8-) ~* * * .

(3.33)

Thus when the Reynolds number is not small, slight coiling reduces the flux. The author (Van Dyke, 1978) has extended this series to ten more terms by computer and by analyzing it deduced that the flux ratio decays for large Dean number as 2 . 1 2 ~ - ’ / ~However, . numerical solutions and several approximate boundary-layer suggest instead a decay as ~ O K - ” ~ . This disagreement remains unresolved. Covering the whole range of Reynolds number requires a double expansion in powers of coiling ratio and Dean number. Larrain and Bonilla have carried that calculation as far as the tenth power. For the first effects on flux ratio they find -Qc =

03

(3.34)

Here the first and third coefficients of E~ appear in (3.31) and (3.33), but the intermediate term is new. Wang (1981) has considered the fact that a coiled pipe is a helix rather than a torus. He includes the first effects of torsion, or pitch, by using a nonorthogonal coordinate system. He finds that the secondary flow is affected significantly at low Reynolds number, but the flux ratio is unchanged to the order of (3.31).

Slow Variations in Continuum Mechanics

35

Earlier, Wang (1980b) used the same coordinate system to calculate the temperature distribution inside an electrically heated circular wire that is loosely coiled into a helix. Perturbing the solution for a straight wire shows that curvature adds corrections of order E = a / L and E’, but torsion again affects only subsequent terms. IV. Three-Dimensional Thin Shapes We remarked in the introduction that a three-dimensional shape with disparate dimensions may be either slender, like a needle, or thin, like a razor blade; and we have so far considered only slender shapes. Thin three-dimensional shapes have been less studied, except for constant thickness as in the theory of elastic plates and shells. We now consider shapes of slowly varying thickness, restricting attention for simplicity to the case of a flat midplane. Then in the notation of Fig. 11 the two faces are described by z = EX, ~ y ) . We saw that for a two-dimensional shape the approximation of slow variations always reduces the mathematical problem to a succession of quadratures. For a slender three-dimensional shape, on the other hand, it is reduced only to partial differential equation in two rather than three variables. Now we shall see that for thin three-dimensional shapes either of these reductions may arise, depending on the nature of the boundary conditions.

FIG. 11.

Symmetric region of slowly varying thickness.

Milton Van Dyke

36

A. STEADYTEMPERATURE I N PLATEOF VARYING THICKNESS The simpler situation-reduction to quadratures-is exemplified by the three-dimensional counterpart of the first problem that we considered in this article: find a solution of Laplace's equation that assumes the values 1 on the two faces of a plate of slowly varying thickness. Again can be interpreted as the steady temperature distribution in the plate when its faces are maintained at different constant temperatures. However, it has now lost its other interpretation as a stream function for irrotational fluid motion, because a three-dimensional flow cannot be described by a single stream function. We now contract y as well as x by introducing X = E X , Y = EY. Then the problem becomes

+

*

+

22

=-

+ +YYL

E"+XX

+

+

=

*1

at

z

=

*f(X, Y).

(4.1)

Quadratures and iteration yield the second approximation $ = ( z / f ) + d."(l/f)xx + (l/f)YYlz(f' - z')

+ O(E4),

(4.2)

and further terms can readily be found. This, of course, reduces (with z replaced by y ) to the two-dimensional result (2.16) when the thickness varies only with x. For comparison, the exact solution for the doubly concave cone having f ( X , Y ) = ( X 2 + Y2)"' is found to be

+ = tanh- ( z / J x 2 + y 2 + z')/tanh-'(E/Jl + Expanding for small E

< 1.

E

E').

(4.3)

shows that in this case the series (4.2)converges for

B. POTENTIALFLOWTHROUGH VARYINGGAP The more complicated situation-reduction to partial differential equations in two variables-is exemplified by potential flow through the slowly varying gap between the two nearly parallel walls of Fig. 11. To a first approximation the velocity components u and ZJ are constant across the gap at any value of x and y, and the transverse component w is negligible. Then the continuity equation becomes

a

-(fu) dX

+ -a( f v ) aY

= 0.

(4.4)

Although we cannot describe the full three-dimensional flow by a single stream function, we can satisfy this quasi-two-dimensional approximation

Slow Variations in Continuum Mechanics

37

by introducing $ according to

fu

fu = +y,

=

(4.5)

-$x.

If the flow is irrotational, all three components of vorticity must vanish, but to a first approximation only the principal component u, - uy normal to the plane of symmetry need vanish. This gives the equation governing $: ($x/fL

+ (Gy/fIy = 0.

(4.6)

This equation is exact for both plane flow (when f = const) and axisymmetric flow (when f = y , the cylindrical radius from the x axis of rotation), but in general it involves a relative error of order E ~ In . polar coordinates, with x = r cos 8 and y = r sin 8, it becomes

As an application we consider uniform potential flow past a sphere tangent to a plane, which is equivalent to flow normal to the line of centers of two equal spheres in contact. This problem was treated by Latta and Hess (1973), who inverted the sphere into a plane and then solved an integral equation involving a Bessel function by a series of transformations. In particular, they evaluated the asymptotic behavior of the flow near the point of contact. Our approximation of slow variation becomes exact as that point is approached, and the spheres can be replaced by the osculating paraboloids. This corresponds to choosing our thickness function asf = X 2 + Y 2= (Er)’, and then (4.7) becomes @rr

-($r/r)

+ ($ee/r2>

= 0.

(4.8)

With the uniform flow directed along 8 = 0, that line can be taken as $ = 0. This suggests trying a solution of the form rC, = r k sin 8, and substituting shows that k = 1 &‘. We choose the less singular solution, which means that $ is locally a multiple of r’+&sin 8. Hence equations (4.5) show that the velocity diverges as rJS-’ = r-0.586approaching the point of contact. This is the result that Latta and Hess found b y more sophisticated means. Sobieczky (1977) has discussed compressible inviscid flow in a slowly varying gap. He considers in particular the transonic small-disturbance approximation for flow through narrow Lava1 nozzles of special form.

*

V. Closer Fits The slow variations discussed so far start from the quasi-cylindrical approximation, according to which the solution at any longitudinal station is that for an infinite cylinder of the local cross section. It is natural to seek

38

Milton Van Dyke

a better approximation by starting with a local solution that fits the boundary more closely. For a boundary of varying width, the next step beyond the quasi-cylinder requires a shape that is at each station tangent to the boundary. For example, in a plane symmetric or axisymmetric problem the solution at each longitudinal station can be approximated by that for the tangent wedge or cone, as indicated in Fig. 12. This approximation requires that only the slope of the boundary be changing slowly rather than its width-in dimensionally correct terms, that the product of width and curvature be everywhere small. In elasticity theory, Massonet (1962) has suggested that the simple solution for the stress field in an infinite wedge loaded at its tip can be used in this way to estimate the distribution of stresses in a beam of variable height. In viscous flow theory, Fraenkel (1962) has selected from the infinite family of self-similar Jeffery-Hamel solutions for laminar flow through a wedge the single one relevant to more realistic channels (the “principal” branch mentioned in our previous discussion of laminar flow in a channel). He has used it (Fraenkel, 1963, 1973) as the basis for an asymptotic theory of flow through a symmetric plane channel with slowly curving walls. His examples show that the higher-order corrections are remarkably small, even when the flow is separated, with regions of reverse motion near the walls. This implies that the first approximation for slowly curving walls is a good one. We observed earlier, in Section I1 (p. 19), that the Jeffery-Hamel solution for a wedge has no axisymmetric counterpart in general, but that self-similar solutions d o exist for an exponential pipe if it is slowly growing. Eagles (1982) has devised a scheme, somewhat analogous to Fraenkel’s, for applying those solutions as a first approximation to a wide variety of “locally exponential” pipes. Again, examples show that the higher-order corrections are very small. However, although reverse flow is included in the basic solutions, it cannot be admitted in the approximation for more general variations of radius.

FIG. 12. Tangent wedge or cone

Slow Variations in Continuum Mechanics

39

The possibility of making an approximation of this sort, starting with a closer fit than the cylindrical one, depends upon the existence of a family of basic solutions, preferably in simple closed form. We have seen that the cylindrical solution can be found for a wide variety of problems. It has a trivial form for the Laplace or biharmonic equation, and for the NavierStokes equations it is known for many cross sections, including ellipses. By contrast, the wedge and cone solutions needed to fit the slope of a boundary are much more limited. For example, not only does the Jeffery-Hamel solution for laminar flow through a wedge have no axisymmetric counterpart for a cone, but even the solution for the wedge involves elliptic functions, which greatly complicates the applications. A. PLANE PROBLEMS

Closer fits are so little developed that in illustrating them we restrict attention to just three plane problems that we have already treated as slow variations. 1. Potential in Symmetric Strip

Consider again the potential in a symmetric strip (Fig. 1). The boundary is no longer restricted to be slowly varying, and so we describe it simply by y = * f ( x ) . (However, a different description might be helpful in higher approximations.) Then according to Fig. 12 the slope f ' ( x ) at any station x is to be identified with the slope E of the wedge in (2.13). Likewise y/(f/y)is to be identified with y / x . Thus the tangent-wedge approximation is found as

+ = tan-'(yf'/f)/tan-'y.

(5.1)

For potential flow through the channel, this gives the velocity along the centerline as u ( x , 0) = +.y(x, 0) = (f'/f>/tan-'f'.

(5.2)

In Table I we compare this approximation with the exact result for a hyperbolic channel of 30" half-opening angle. We also show for comparison the first approximation of slow variations, from (2.5). We might have anticipated that the closer fit would always prove more accurate. However, the result of slow variations is seen to be the more accurate near the throat of the channel, where the slope of the boundary is small, the tangent-wedge approximation being superior only for x > 2, where the product of wall curvature and channel width is small. We saw in Fig. 3 how the solution of this problem is improved by extending the method of slow variations to higher order. It would be desirable in the

Milton Van Dyke

40 VELOCITY

AXIS

ON

X

Slope

0 0.5 1.o 1.5 2 3 4 6 10

0 0.160 0.289 0.378 0.436 0.500 0.530 0.555 0.569

TABLE I POTENTIAL FLOW THROUGH HYPERBOLIC CHANNEL y = + ( I + x2/3)"' FOR

Curvature x width

Tangentwedge

Slow variation

Exact

0.666 0.569 0.384 0.235 0.144 0.060 0.029 0.010 0.002

1.000 0.969 0.890 0.791 0.694 0.539 0.432 0.304 0.188

1.000 0.961 0.866 0.756 0.655 0.500 0.397 0.277 0.171

0.955 0.926 0.854 0.764 0.675 0.530 0.427 0.302 0.187

same way to embed the tangent-wedge method into a systematic scheme of successive approximations. However, it is not clear how this is to be done. For plane laminar flow in a symmetric channel Fraenkel (1963), starting with the Jeffery-Hamel solution for flow in a wedge, has carried out the third approximation. However, he requires that the channel be mapped conformally onto a strip, which in our present example of potential flow means knowing the exact solution. 2. Potential in Meandering Strip

We now reconsider the meandering strip of constant width shown in Fig. 6 b . Whereas we previously assumed that the curvature K of the centerline was small as well as slowly varying by setting K = E g ( E s ) , we now assume it to be slowly varying but not necessarily small, with K

=

G(Es).

(5.3)

Then Laplace's equation ( 2 . 5 3 ) assumes, instead of (2.55), the only slightly more complicated form

a

-(1+ an

With terms in

E*

a+

Gn)-+ an

a 1 a4 = 0. aSl+GnaS

E*---

(5.4)

neglected, quadratures yield the first approximation +=1-2

ln[(l

+ G)/(1 + G n ) ]

ln[(l + G ) / ( 1- GI1

(5.5)

Comparison with (2.59) shows that this is the exact solution for a circular annulus that at each station fits the local curvature of the strip, as indicated in Fig. 13. We may call it the osculating annulus approximation.

Slow Variations in Continuum Mechanics

41

FIG. 13. Osculating annulus approximation for meandering strip.

In this case it is clear how higher approximations can be found. Substituting our first approximation (5.5) into the neglected terms in (5.4) yields an iteration equation that can again be solved by quadratures. The resulting correction of order E’ is found in terms of elementary functions-powers of (1 + G n ) and its logarithm (and this would appear to be true of all subsequent approximations). However, the results is too complicated to be given here. 3. Laminar Flow in a Meandering Channel

This osculating annulus approximation can be applied also to plane laminar flow (Van Dyke, 1983) because the required basic solution exists in simple closed form (2.60). (By contrast, we are unable to treat in this way the three-dimensional flow through a meandering pipe of circular cross section, because except at zero Reynolds number the required basic solution for laminar flow in a torus can be found only approximately, as shown in Section III,B,2.) a. Reynolds Number of Order Unity If the Reynolds number is fixed, introducing the slowly varying centerline curvature G ( E sof ) (5.3) into the vorticity equation (2.60) shows that the first approximation satisfies

a

-(1+

an

a i a alC, Gn)--(1 + Gn)- = 0. an dn 1 + Gn an

Milton Van Dyke

42 Four quadratures yield

+ = A( 1 + G n ) 2In( 1 + G n ) + B(1 + G n ) 2+ C In( 1 + G n ) + D.

(5.7)

this is just our previous solution (2.62) for laminar flow in an annulus except that now the coefficients A, B, C, D, which are determined by imposing the boundary conditions at the walls, depend on the axial coordinate through the slowly varying curvature function G (E S ) . The neglected terms are of order E R and E ~ and , higher approximations would proceed in powers of those two quantities. Of course, the analysis is too complicated to be carried much farther by hand, but the author has suggested that it could be continued by using a computer program that manipulates symbols, such as MACSYMA. b. High Reynolds Number When E R is of order unity, the first approximation is governed, with an error of order E’, by the nonlinear equation [;(l+Cn)--ER an a

-----

(::a:

::a:)]

a a* (1 + G n ) 1 + Gn a n an i

~-

= 0. ( 5 . 8 )

This represents those terms in the Navier-Stokes equations that contribute to second-order boundary-layer theory, with its effects of longitudinal curvature. This equation, like its counterpart (2.63) for the classical boundary layer, is parabolic, with no upstream influence, and can therefore be integrated numerically along the channel (at least so long as no reverse flow occurs). Blottner (1977) has carried out the calculations for a diverging channel with a semicircular centerline. When the divergence is enough to produce a region of reverse flow at the inner wall, the numerical integration becomes unstable shortly beyond the separation point. VI. Concluding Remarks

The approximation of slow variations has not received the attention that it deserves. Our survey of the method shows that it has been only sporadically developed, except perhaps for steady laminar flows. Further exploitation will surely yield many informative and useful results in the various branches of continuum mechanics. Higher approximations will be needed, both to refine the accuracy and to establish limits of convergence. Hand calculation becomes infeasible at an early stage, even in simple problems; but the routine labor involved can be readily delegated to a computer. Lucas’s (1972) treatment of laminar flow through a channel shows how even numerical computation can be

Slow Variations in Continuum Mechanics

43

used effectively. Further advances are to be anticipated from the application of symbol-manipulation programs. No doubt the majority of useful results will be based on the simple quasi-cylindrical approximation, because that can be found in simple form for a wide variety of problems in many fields. However, serious consideration must be given to the suggestion of Fraenkel (1962) for laminar flow that greater accuracy will result when a basic solution can be found that fits the geometry more closely. The technique for calculating higher approximations on that basis will have to be developed, for it was not apparent to us how to improve systematically even the tangent-wedge approximation for plane potential flow through a symmetric channel. Similarly, a systematic technique must be developed for calculating higher approximations in thin three-dimensional regions. We saw that for potential flow through the gap between two nearly parallel walls a quasi-twodimensional stream function provides a convenient first approximation, but it is not obvious how that is to be refined. A lot of good hard thinking needs to be devoted to the idea of slow variations ! ACKNOWLEDGMENTS The writing of this article was supported by the National Science Foundation under Grant No. MSM-8300537. The author is indebted to Andreas Acrivos for helpful discussion. REFERENCES Abramowitz, M. (1949). On backflow of a viscous fluid in a diverging channel. 1.Marh. Phys. (Cambridge, Mass.) 28, 1-21. Blasius, H. (1910). Laminare Stromung in Kanalen wechselnder Breite. 2. Marh. Phys. 58, 225-233. Blottner, F.G. (1977). Numerical solution of slender channel laminar flows. Comp. Methods Appl. Mech. Eng. 11. 319-339. Chow, J. C . F., and Soda, K. (1972). Laminar flow in tubes with constriction. Phys. Fluids 15, 1700- 1706. Daniels, P. G., and Eagles, P. M. (1979). High Reynolds number flows in exponential tubes of slow variation. J. Fluid Mech. 90, 305-314. Dean, W. R. (1928). The stream-line motion of fluid in a curved pipe. Philos. Mag. 5(7), 673-695. Eagles, P. M. (1982). Steady flow in locally exponential tubes. Proc. R. Soc. London Ser. A 383, 231-245. Eagles, P. M., and Muwezwa, M. E. (1986). Approximations to flow in slender tubes. J. Eng. Math. 20, 51-61. Eagles, P. M., and Smith, F. T. (1980). The influence of nonparallelism in channel flow stability. J. Eng. Math. 14, 219-237. ErdClyi, A. (1961). An expansion procedure for singular perturbations. Atti Accad. Sci. Torino Cl. Sci. Fis. Mat. Nat. 95, 651-672. Fraenkel, L. E. (1962). Laminar flow in symmetrical channels with slightly curved walls. I. On the Jeffery-Hamel solutions for flow between plane walls. Proc. R. SOC.London Ser. A 267. 119-138.

.(

44

Milton Van Dyke

Fraenkel, L. E. (1963). Laminar flow in symmetrical channels with slightly curved walls. 11. An asymptotic series for the stream function. Proc. R. SOC.London Ser. A 272, 406-428. Fraenkel, L. E. (1973). On a theory of laminar flow in channels of a certain class. Proc. Cambridge Philos. SOC.73, 361-390. Geer, J., and Keller, J. B. (1979). Slender streams. J. Fluid Mech. 93, 97-115. Gohner, 0. (1930). Schubspannungsverteilung im Querschnitt einer Schraubenfeder. Ing. Arch. 1, 619-644. Gohner, 0. (1931a). Schubspannungsverteilung im Querschnitt eines gedrillten Ringstabs mit Anwendung auf Schraubenfedern. Ing. Arch. 2, 1-19. Gohner, 0. (1931b). Spannungsverteilung in einem an den Endquerschnitten belasteten Ringstabsektor. Ing. Arch. 2, 381-414. Kaimal, M. R. (1979). Low Reynolds number flow of a dilute suspension in slowly varying tubes. Inf. J. Eng. Sci. 17, 615-624. Keller, J. B., and Geer, J. F. (1973). Flows of thin streams with free boundaries. J. FZuid Mech. 59, 417-432. Kotorynski, W. P. (1979). Slowly varying channel flows in three dimensions. J. Insf. Mafh. Its Appl. 24, 71-80. Larrain, J., and Bonilla, C. F. (1970). Theoretical analysis of pressure drop in the laminar flow of fluid in a coiled pipe. Trans. SOC.Rheol. 14, 135-147. Latta, G. E., and Hess, G. B. (1973). Potential flow past a sphere tangent to a plane. Phys. Fluids 16, 974-976. Lucas, R. D. (1972). A perturbation solution for viscous incompressible flow in channels. Ph.D. dissertation, Stanford Univ.; Univ. Microfilms, order no. AAD72-30664. Manton, M. J. (1971). Low Reynolds number flow in slowly varying axisymmetric tubes. J. Fluid Mech. 49, 451-459. Massonet, Ch. (1962). Elasticity: Two-dimensional problems, In “Handbook of Engineering Mechanics” (W. Fliigge, ed.), pp. 37-1 to 37-30, especially p. 37-23. McGraw-Hill, New York. Morse, P. M., and Feshbach, H. (1953). “Methods of Theoretical Physics.” McGraw-Hill, New York. Olson, D. E. (1971). Fluid mechanics relevant to respiration: Flow within curved or elliptical tubes and bifurcating systems. Ph.D. thesis, Imperial College, London. Rosenhead, L., ed. (1963). “Laminar Boundary Layers.” Oxford Univ. Press, London. Sobey, I. J. (1976). Inviscid secondary motions in a tube of slowly varying ellipticity. J. Fluid Mech. 13, 621-639. Sobieczky, H. (1977). Kompressible Stromung in einer ebenen Schicht variabler Dicke. Z. Angew. Mafh. Mech. 57, T 207-T 209. Timoshenko, S., and Goodier, J. N. (1951). “Theory of Elasticity.” McGraw-Hill, New York. Todd, L. (1977). Some comments on steady, laminar flow through twisted pipes. J. Eng. Math. 11, 29-48. Todd, L. (1978). Steady, laminar flow through non-uniform, curved pipes of small cross-section. Tech. Rept. No. 78-19, Inst. Appl. Math. Stat., Univ. British Columbia, Vancouver. Todd, L. (1979). Steady, laminar flow through curved pipes of small, constant cross-section. Tech. Rept. No. 79-1, Dept. Math., Laurentian Univ., Sudbury, Canada. Todd, L. (1980). Steady, laminar flow through non-uniform, thin pipes. Tech. Rept. No. 80-1, Dept. Math., Laurentian Univ., Sudbury, Canada. Van Dyke, M. (1970). Entry flow in a channel. J. Fluid Mech. 44, 813-823. Van Dyke, M. (1975). “Perturbation Methods in Fluid Mechanics,” Parabolic, Stanford, California. Van Dyke, M. (1978). Extended Stokes series: Laminar flow through a loosely coiled pipe. J. Fluid Mech. 86, 129-145. Van Dyke, M. (1983). Laminar flow in a meandering channel. SIAM J. Appl. Math. 43,696-702.

Slow Variations in Continuum Mechanics

45

Wang, C.-Y. (1980a). Flow in narrow curved channels. J. Appl. Mech. 47, 7-10. Wang, C.-Y. (1980b). The helical coordinate system and the temperature inside a helical coil. J. Appl. Mech. 47, 951-953. Wang, C. Y. (1981). On the low-Reynolds-number Row in a helical pipe. J. Fluid Mech. 108, 185-194. Wild, R., Pedley, T. J., and Riley, D. S. (1977). Viscous Row in collapsible tubes of slowly varying elliptical cross-section. J. Fluid Mech. 81, 273-294. Williams, J . C., 111 (1963). Viscous compressible and incompressible Row in slender channels. A I M J. 1, 186-195.