LONGITUDINAL DISPERSION FOR TURBULENT
FLOW IN PIPES C. N. SITTEL, JR., W. Vanderbilt University, Nashville, Tenn.
K . B. S C H N E L L E , J R .
D. THREADGILL, AND
37203
The applicability of the one-dimensional turbulent dispersion equation as a model to describe longitudinal dispersion in pipe flow was investigated by pulse testing for the NaCI-water system in test sections of 72 feet of glass pipe of different diameters and with conductivity measuring cells located at varying distances from the base cell. The model was valid for describing longitudinal dispersion for Reynolds numbers greater ihan 4 X 1 04. In fully developed turbulent flow the dispersion coefficient was independent of pipe length; the effects of diameter and fluid velocity were in agreement with the model.
turbulent disperr,ion model has been used in studies of tubular reactors, packed-bed extraction columns, chromatography columns, heat exchangers, and pollutant loads in natural streams. T h e objectives of this study were to investigate the applicability of the one-dimensional turbulent dispersion equation to describe longitudinal dispersion in pipe flow and the effect of pipe diameter, length, and fluid velocity on the dispersion coefficient; and to compare the dispersion coefficients obtained from this model with those predicted by Taylor (19 54). HE
Literature Review
Solutions of the partial differential equations describing turbulent diffusion in one dimension have been presented for a variety of boundary conditions by Bischoff and Levenspiel (1962, 1964), Clements (1963), Croockewit et al. (1955), Hays (1964), Levenspiel (1904), Levenspiel and Smith (1957), and Wehner and Wilhelm (1956). Taylor (1954) derived an expression for the virtual coefficient of diffusion
K = 10.1
TU
little difference between the dispersion coefficients in 12- and 10-inch diameter systems. Theory
In a turbulent flow system through which a tracer is carried by a homogeneous solvent, bulk flow and longitudinal dispersion affect the distribution of the tracer in space and time. If the effects of radial dispersion, velocity profile, and mass transport by molecular diffusion are assumed negligible and the velocity and dispersion coefficient are considered to be constant, a differential mass balance gives
which is the basic equation for turbulent dispersion that was investigated in this study. Substituting the dimensionless variables 21
e = - t L
)(:
z
= x/L
Pe = uL/DL
where in Equation 2 gives by considering the effects of radial diffusion, velocity profile, and the analogy between momentum and mass transfer in turbulent flow. Levenspiel and Bischoff (1962, 1964) summarized the mathematical models used in describing flow in a number of physical situations and presented methods for determining the Peclet number and the mean residence time from the mean and variance of the pulse tracer curves. For axial dispersion, they found that the data were represented by a straight line on a logarithmic plot of the ratio of the dispersion coefficient to the kinematic viscosity us. the Reynolds number. Davidson et a l . (1955) and Hull and Kent (1952) measured the dispersion of petroleum products flowing through long pipelines. These studies were not performed under desirable conditions, since Taylor (1954) had noted that pipe bends and fittings altered tbe magnitude of the dispersion coefficient. According to Davidson et a i . (1955), the Schmidt number had no effect on the dispersion coefficient. Schlinger and Sage (1953), using diesel fuel and naphtha as tracer with air as carrier fluid, fqund that a t equal Reynolds numbers there was
(3) Clements’ (1963) solution of this equation for input and output traces measured inside the system was used in this work. T h e boundary conditions for this solution were: C(z,O) = 0
Lim [C(z,e)] = 0 z+m
C(0,e) = C&) for e
>0
T h e points of measurement of the input and output pulses were z = 0 and z = 1, respectively. T h e Laplace transform of Equation 3 with respect to dimensionless time is
This solution was divided into two parts: VOL. 7
NO. 1
FEBRUARY 1968
39
CI(S),the Laplace transform of the input pulse and
f ( z , S ) , the Laplace transform of the system function, giving
ments of the pulse, composed of three equally spaced time points, were fitted by a parabolic equation. T h e first two segments may be represented mathematically by the quadratic expressions
yi(t) = ait2 f bit f T h e inversion of this solution into the time domain would require numerical convolution of the transformed functions. T h e difficulty involved in numerical convolution, however, could be reduced by using the convolution theorem (Sokolnikoff and Redheffer, 1958), transforming the two functions into the frequency domain, and then multiplying them together. T h e Fourier transform for the system function was obtained by substituting the dimensionless frequency, iW, for S in Equation 4. Thus, F(z,W) = exp
[I
- VI'
(6)
Defining
4 W/Pe
and b = l/z tan-' a and applying de Moivre's theorem, the expression for F(z, W) becomes F(z,W) = exp
{F
- 11
0
5 t < Ti
(14)
respectively. T h e Laplace transforms for the second derivatives ofyl(t) and y z ( t ) are, respectively,
and
f 4 iw/pe] 2/Pe
a =
CI
and
The first and second bracketed terms represent additions to the transform caused by discontinuities in the slope and in the second derivative, respectively. T h e exponential term gives the time elapsed between zero and the start of the particular segment. I n transforming into the frequency domain, iw was substituted for s and the exponential expanded into complex form by Euler's relation. T h e real, R, and imaginary, Z, parts for 1)th approximating equations were: the (i
+
[.(+I
- (1 f a2)1~4cos(b)i(l f az)l/4sin(b)]}
- a t ] sin (Tim) -
(7) and
This equation rewritten in terms of magnitude ratio, M.R., and phase angle, P.A., is
=
w2
F(z,W) = [M.R.][cos (P.A.) f i sin (P.A.)]
-
{2w
[U{+l
- a * ] cos ( T @ )-
(8)
where T h e real and imaginary parts for the entire pulse were obtained by summing over all approximating sections for each frequency. A linear relation of the form
and P.A. =
- Pe - (1 9a2)1/4sin(6) 2
c = mgfp
For computational purposes, the system function was converted into real, RM,and imaginary, ZM, parts:
RM
=
M . R . COS (P.A.)
(11)
ZM
=
M.R. sin (P.A.)
(12)
and
I n calculating the Fourier transforms of the input and output pulses, the method described by Hays et al. (1965) was used. T h e general relation for a Fourier transform is
T [ g ( t ) ]=
g ( t ) e-"'''
Jm -m
dt
(1 3)
Since the pulse is nonexistent prior to time zero, the integral with limits from minus infinity to zero is zero. Aseltine (1958) noted that for this situation, the Laplace transform was equal to the Fourier transform with s substituted for iw. T o derive a numerical transform of the tracer curve, seg40
IhEC FUNDAMENTALS
(20)
existed between concentration and conductivity. Since the concentration term in Equation 3 is due solely to the tracer and since the water contained dissolved impurities, a fictitious base tracer concentration, c,, corresponding to the equivalent concentration of the impurities was subtracted from the measured concentration, giving
or
Each cell had a different cell constant and it was necessary to normalize the data in processing. T h e mass of tracer in the system, represented by the area under the pulse curve, was assumed to be constant. Thus the normalized conductance was defined as
(23)
Equation 3 rewritten in terms of the normalized conductance is
bG
be
-
1 dZG Pe dz2
_ -bG
(24)
dz
\vith boundary conditions G(r,O:i = 0 Lim [G(z,e)]= 0 z+a
G(0,e:r
= G,(e)
for e
>0
Therefore the only experimental data needed for testing the mathematical model were conductivity measurements. I n a fitting procedure for obtaining the statistically best values of a parameter, the integral of the squared error is the function which is usually minimized. This function is defined as
This equation was examined for errors introduced in transforming into the frequency domain. If the time domain error, e ( t ) , is defined as
4 ) = .ro(t) - Y p ( Q
(26)
the expression for 9 mav be written as 9
==
[ e ( t )I2dt
Jm
(27)
Parseval's theorem (Lepage. 1961) applied to Equation 27 gives
Lrn
(28)
where the Fourier transform of e ( t ) is E&).
E(iw) is defined as
9 = J m [e(t)lz dt =
E(tw) == Y,(iu)
A
'E(iw)
- Yp(iu)
Experimental Equipment
(29)
where Y,(iw) and Y p ( i u )are the frequency transforms of Yo(t) and Y p ( t ) respectively. , There is a direct correlation between a minimization of the squared error and a minimization of the absolute squared transfcrmed error in the frequency domain. T h e predicted output consists of a complex multiplication of the transformed input, R , and I,, by the frequency expression of the model, RII and ZAPI; thus
Y,(iw) = (R,
+ 2,) ( R w + i l w )
(30)
T h e complex expression for the observed output pulse is
Y,(zu)
=
R,
+ iZ,
Sololnikoff and Redheffer (1958) state that the partial sums of a Fourier series give a smaller mean-square error than any other available linear combination, and that successive terms in the Fourier series tend to zero as the number of terms increases. These statements imply that the model has been approximated by the best least-squares method available, and that if w c has been chosen sufficiently large the loss in accuracy is small. Truncation of the frequency spectrum also produces a beneficial effect, in that it discards the erroneous fluctuations in the transformed data a t high frequencies, caused by inaccuracies in reading the data and by use of parabolic segments to represent the continuous curves. Equation 33 was minimized using the nonlinear least squares program of Marquardt and Baumeister (1962) and the statistically best values of the Peclet number and the mean residence time, L / u . Bischoff and Levenspiel (1962) showed that the first moment of the tracer curve about the origin provides a good estimation of the average residence time, L/u. T h e Peclet number was estimated from a study of the local minima surrounding the true minimum. These two estimates were used as the initial values in the nonlinear least-squares program. In the numerical regression analysis, the problem of local minima existed. T o determine their location the sum of the squared error, 9, was plotted against the Peclet number for each iteration from the nonlinear least squares fitting of the first cell against the second one. When the Peclet number was chosen below or slightly above the true minimum, the solution proceeded to this true minimum, but when the estimate was considerably above the true minimum, the results were those from local minima. T o ensure that all the results were those from the true minimum. the sum of the squared error, 9,was plotted against the Peclet number for the largest and one of the smallest flow rates for each pipe diameter investigated. Values of Peclet number for intermediate flow rates were obtained by interpolation or extrapolation.
T h e apparatus (Sittel, 1966) was similar to that used by Hays (1964). Conductivity-measuring cells were placed 10, 20, 30, 50, and 70 feet from the input concentration-measuring cell. The test section was 72 feet long, constructed from I/*-, 1-, and 2-inch nominal diameter glass pipe. Ten feet of black iron pipe were placed before the test section for calming, and the first cell was located 18 inches downstream from the entrance to the test section. One 6-inch section of glass pipe was placed between the last cell and the drain line. T a p water was pumped through the apparatus by a centrifugal pump. T h e test apparatus and the injection system are shown schematically in Figure 1. I t was desirable to introduce the tracer as rapidly as possible to ensure a high frequency content of the input pulse (Clements
(31)
COMPRESSED
T h e complete expression for 9 is 1
TRACER
RESERVOIR
/.-
SOLENOID VALVE
CLASS PIPE
1 / / /II
In processing the data om a digital computer, the relation for CP was made discrete. L,etting Z ( u n ) represent the integrand in Equation 32, the relation for 9 becomes l
m
A
n=t
Z(w,)
9 ==
I kI0'
+ I
I
I
I I/
10'+-20'+20'+
DISTANCE BETWEEN CELLS
DRAIN
CITY WATER MAIN
(33)
This summation was truncated a t a finite frequency limit, w c .
CENTRIFUOAL
PUMP
Figure 1.
Schematic of test apparatus
VOL. 7
NO.
1
FEBRUARY 1968
41
and Schnelle, 1963). T o ensure cross-sectional homogeneity of the tracer solution a t the first cell, the injection port was placed upstream from this cell. Although spreading of the tracer in this section reduced the frequency content of the measured input pulse, cross-sectional homogeneity was necessary, since the mathematical model was derived for dispersion in only one dimension. The injection apparatus was a solenoid valve with a 1/8-inch exit port drilled through the base. A Teflon gasket, coated with cement, ensured water-tight seal between the base of the valve and a 1/8-inch injection port bored in the calming section. A 5-gallon pressure vessel fed the tracer to the solenoid valve. T h e contents of the tank were maintained a t 140 p.s.i.g. An Agastat time-delay relay actuated the solenoid valve for a satisfactory duration of injection. I n turbulent flow the effect of molecular diffusion was negligible compared to longitudinal dispersion (Davidson et al., 1955; Taylor, 1954). Therefore, any tracer which had a linear relationship between concentration and conductivity could be used. Solutions of sodium chloride in tap water were tested, and the linearity was such that additional calculations were not needed to process the data by Equation 5. A reduction of the base resistance by a factor of Z1/2 to 4 by the tracer was desired, because this represented a region of linearity between concentration and conductivity and was within the limits of the apparatus. Saturated stock solutions of table salt in water diluted to one part of stock to 25 parts of water for the 11/2-inch diameter test section and 1 to 1 5 for the 1- and 2-inch diameter test sections gave the desired decrease in resistance. T h e conductivity cells were constructed from Teflon-T gaskets and plated with platinum black to reduce polarization. Holes were drilled through the edges so that a maximum crosssectional area was available for measuring conductivity. Twenty-gage platinum wire was inserted through the holes, and the entire cell was cleaned by immersion in boiling aqua regia. A I-kc. sine wave was fed to a power amplifier with source impedance of 8 ohms. This voltage source was connected in series to each cell and its corresponding measuring bridge. T h e voltage was rectified for recording purposes by a full wave rectifier. An oscillographic recorder was used. Results
The relationship between dispersion coefficient and Reynolds number obtained from nonlinear least squares fitting of the 387 data points is
D L = 3.87 X 10-5 Re0.764, sq. ft./sec.
(34)
This equation and the experimental data are shown graphically in Figure 2, where each point represents an average value per run obtained from the values for each of the five lengths used. For comparison the results of Hays (1964) and Taylor (1954) are included. Calculation of u* by combination of the Blasius relation f=-
0.0791 Re1i4
=
I
t
1
I
I
I
1
I
1
Ob
3 0.3
< 0.2 W
N.
- 0.1 U
n ’
.06
Figure 2. Effect of Reynolds number on average dispersion coefficient
fore, it was concluded that the results of this study were not significantly different from those of Taylor. Figure 3 indicates the ability of the model to fit the data and the total error involved. T h e standard error of estimate,
was the criterion used in determining the “goodness of fit” of the model. Figure 3 indicates that the standard error of estimate tends to increase as the Reynolds number decreases. This is as expected, since the model was derived for fully developed turbulent flow. At Reynolds numbers greater than 4 X lo4 the standard error of estimate was approximately constant. For fully developed turbulent flow the standard error of estimate should be zero; therefore the deviation from zero is a measure of the experimental and model errors. For Reynolds numbers less than 4 X lo4 the major contribution to the standard error of estimate was probably the fact that the turbulent flow assumption in the model was invalid. The results indicate that the dispersion coefficient is essentially independent of length except for low Reynolds numbers (less than 4 X lo4). Since the Model was derived for fully developed turbulent flow, the coefficients calculated a t low Reynolds numbers would be expected to be in error. Runs were made to determine the effect of the measuring cells on the dispersion coefficient. The Reynolds number was held as constant as possible and the system was pulsed with all cells in place, with the second cell removed, with the second and third cells removed, and finally with only the first and last cells in place. If the presence of the platinum wires in the flow pattern significantly affected the turbulence, the dispersion
(35)
o’
and the relation (Bird et al., 1960)
f
I .o
0.6
2(:)z
and substitution in Taylor’s equation (Equation 1) give for this system
D L = 1.22 X 10-5 Re0.875, sq. ft./sec.
(37)
Ninety-nine per cent confidence limits for Equation 34 sq. foot per second for the are -2.82 X 10-5 to 1.06 X coefficient term, and 0.618 to 0.910 for the power term. The constants in Taylor’s equation are within these limits; there42
I L E C FUNDAMENTALS
0.2
0.oL
IO3
I
1
1
2 3
6
IO4
I
I
I
I
2
3
6
1
I
I
I
1
IO5
2
3
6
10‘
Re Figure 3.
Effect of Reynolds number on total error
=
coefficient should decrease as cells are removed because of the decrease in turbulence. No significant change in dispersion coefficient was found. To determine whether changes in flow rate affected the conductance measurements, the system was purged and the floiv rate was varied over the attainable range. No change in base conductance was detected.
= = = =
= = = =
Conclusions
=
T h e one-dimensional turbulent dispersion equation is valid for describing longitudinal dispersion for Reynolds numbers greater than 4 X lo4. In fully developed turbulent flow, the longitudinal dispersion coefficient is independent of length. T h e effects of pipe diameter and flow rate on the longitudinal dispersion coefficient are in agreement with the onedimensional equation. Acknowledgment
a,
c c,
= = = = =
CO
=
b b,
C
Cr(S) CI(0)
= =
=
C(z,S) =
DL
=
e(t)
=
E(iw)
= =
f(z,S‘)
= = =
f
F(t,H’) g go g(t)
G I
= = = =
I.,f
= =
IO
=
I
=
K L
=
I1
m M.R. n
P
P.A. Pe r
=
= = = = = = =
= = = = = =
= = = = = = = =
=
Nomenclature = 4 11’/Pe
= =
=
C. E.Sittel, Jr.. gratefully acknowledges financial support in the form of a traineeship from the National Aeronautics and Space Administration.
a
v* W
=
=
polynomial coefficient 1/2 tan+ a polynomial coefficient concentration of tracer, lb. moles/cu. ft. polynomial coefficient base concentration. lb. moles/cu. ft. effective concentration of tracer, lb. moles,/cu. ft. Laplace transform of input pulse, lb. moles/cu. ft. concentration of input pulse in dimensionless time, lb. moles/cu. ft. Laplace transform solution of model, lb. moles/cu. ft . longitudinal or axial dispersion coefficient, sq. ft./ sec. error term Fourier transrorm of error term Fanning friction factor Laplace transform of system function Fourier transform of system function conductance. mhos base conductance, mhos arbitrary function in time normalized conductance, set.-' imaginary part for i l t h approximating equation of pulse imaginary part of input pulse imaginary part of mathematical model imaginary part of output pulse denotes imaginary root, 1/-1 virtual coefficient of diffusion distance between initial and final points of measurement, ft. polynomial coefficient magnitude ratio summation index polynomial coefficient phase angle Peclet number, u L / D L pipe radius
+
= =
+
real part of i l t h approximating equation of pulse real part of input pulse real part of mathematical model real part of output pulse Reynolds number, d u p / p Laplace transform variable dimensionless Laplace transform variable time, sec. time a t end of ith approximating section, sec. velocity, ft./sec. turbulent velocity, ft./sec. dimensionless frequency, Lw/u distance variable, ft. Laplace transform of first approximating equation first approximating equation Laplace transform of second approximating equation second approximating equation observed output predicted output Fourier transform of observed output Fourier transform of predicted output dimensionless length, x / L squared error in frequency domain effective conductance, mhos dimensionless time, ut/L sum or integral of squared error wall shear stress lb./ft. s e c 2 frequency, sec-1 cutoff frequency, sec.-1 frequencies over which summation is made, sec.-l denotes convolution density, lb./cu. ft.
literature Cited
Aseltine, J. A., “Transform Method in Linear System Analysis,” McGraw-Hill, New York, 1958. Bird, R. B., Stewart, W.E., Lightfoot, E. N., “Transport Phenomena,’’ LViley, New York, 1960. Bischoff, K. B., Levenspiel, O., Chem. Eng. Scz. 17, 245, 257 (1962). Clements, i V . C., Jr., “Pulse Testing for Dynamic Analysis,” Ph.D. thesis, Vanderbilt University, 1963. Clements, \$’. C., Jr., Schnelle, K. B., Jr., Znd. Eng. Chem. Process Design Decelop. 2, 94 (1963). Croockewit, P., Honig, C. C., Kramers, H., Chem. Eng.Sci. 4, 111 (1955). Davidson, F. F., Farquharson, D. C., Picken, J. O., Taylor, D. C., Chem. Eng. Scz. 4,201 (1955). Hays, J. R., “Mathematical Modeling of Turbulent Diffusion Phenomena,” M.S. thesis, Vanderbilt University, 1964. Haw, J. R., Clements, LV. C., Jr., Harris, T. R., “Model Evaluation in the Frequency Domain,’’ 58th Annual Meeting, A.I. Ch.E., Symposium on Mathematical Modeling of Chemical Processes, December 1965. Hull, D. E., Kent, J. FV., Ind. Eng. Chem. 44,2745 (1952). Lepage, \V. R., “Complex Variables and the Laplace Transform for Engineers,” McGraw-Hill, New York, 1961, Levenspiel, O., “Chemical Reaction Engineering,” Wiley, New York, 1964. Levenspiel, O., Bischoff, K. B., Aduan. Chem. Eng. 4, 95 (1964). Levenspiel, O., Smith, W.K., Chem. Eng. Sci. 6,227 (1957). Marquardt, D. LV., Baumeister, T., Share General Program Library, Share Distribution No. 1428 DFE 2135 (PA) (1962). Schlinger, W.G., Sage, B. G., Znd. Eng. Chem. 45,657 (1953). Sittel, C . N., Jr., “Longitudinal Dispersion for Turbulent Flow in Pipes,” M.S. thesis, Vanderbilt University, 1966. Sokolnikoff, I. S., Redheffer, R. N., “Mathematics of Physics and Modern Engineering,” McGraw-Hill, 1958. Taylor, G . I., PTOC. Roy. SOC. (London) A223, 446 (1954). LVehner, J. F., Wilhelm, R. H., Chem. Eng. Sci. 6,89 (1956). RECEIVED for review January 3, 1967 ACCEPTED August 7, 1967
VOL. 7
NO. 1
FEBRUARY 1968
43
