- Email: [email protected]
A puff pollutant dispersion model with wind shear and dynamic plume rise
Atmos+rk Entdmmmt Vol. If pp. 1933-1938. 0 Pwpttton Press Ltd. 1978. Prinral in Gmt FJriuin.
A PUFF POLLUTANT DISPERSION MODEL WITH WIND SHEAR AND DYNAMIC PLUME RISE* C. M. Atm~phe~c
SHEIH
Physics Section, Radiolo~~l and ~vironm~tal Research Division, Argonne National Laboratory, Argonne, IL 60439, U.S.A.
(First receiwd 8 Nuuember 1977 and in&al form 6 March 1978)
Abstract - A puff diffusion model, which includes wind shear and dynamic plume rise, is developed for numerical prediction of pollutant concentrations under unsteady and non-uniform flow conditions. The plume from a continuous source is treated as a series of put%emitted successively from the source. Each p&is represented by a set of six tracer particles, which define the size, shape and location of the puff. Initially these. particles are located at the surface of the source, on arbitrarily chosen orthogonal axes. The location of the particles is computed at each time step by taking into account advection, eddy diffusion, wind shear and entrainment of ambient air during plume rise. The concentration distribution of each puff isdetermined by fitting an ellipsoid to thecluster of the six particles and assuming a three-dimensional Gaussian distribution, with standard deviation quaf to the half-lengths of the principaf axes of the ellipsoid. The concentration at a point of interest is obtained by summing the ~n~butions from nearby puffs. The effect of wind shear on the pollutant con~tration is investigated by use of a typical wind shear encountered in the atmosphere. The resufts show that, at 6OOmdownstream from the source, the present mode1gives con~trations a factor of 2 higher and lower at one standard deviation below and above the plume center, respectiveiy, than that of conventional models in which no wind shear is considered. The phtme&e formulation is calibrated against the observations compiled by Briggs and the model is used to predict the trajectory of a plume observed by Slawson and Csanady. Excellent agreement between the prediction and the observation can be achieved if an appropriate eddy diffusivity is chosen.
1. iNTRODUCTION
A considerable number of diffusion models have been developed during the past several decades. Basically, the available models can be classified into plume, puiI and grid models. Plume models have been developed and applied by many investigators, such as Lucas (1958), Turner (19&Q, Hilst (1%7), and Shieh er al. (1970). In these models the growth of the plume width is a function of atmospheric stability and of travel time, as summarized by Turner (1969), and the concentration distribution of the plume in the cross-wind direction is assumed to be Gaussian. However, plume models are derived under the assumption of spatial and temporal uniformity in wind and atmospheric stability. Obviously, these assumptions are rather limiting for many air pollution problems in which factors such as non-~ifo~ to~graphical features, sea or lake breeze circulations, and urban heat-island effects can cause significant temporal and spatial variations in meteorological conditions. In order to include local effects, Lamb (1%9), Roberts er al. (1970) and Ludwig et al. (1977) have developed puff models. In these models, source emissions are treated as a series of puffs emitted into the atmosphere The concentration distribution within each pnfI is assumed to be Gaussian. In neither plume nor puff models is wind shear taken into account. This can be a rather serious * Work supported by U.S. Department of Energy and The Pennsylvania State University.
omission, as the observations by Schiermeier and Niemeyer (1970) show that plume cross-sections are highly elongated by wind shear. Although the e&t of wind shear on an isolated puff has been investigated by Richards (197(l), the results can not be applied directly to the present problem because the present puff is not an isolated one but a segment of a continuous plume. Moreover, Richards’ solution is obtained under constant shear and constant total buoyancy force, which is inappropriate for the non-uniform, unsteady state addressed in the present paper. In general, grid models solve finite-difference equations of diffusion. These models have the potential to be more realistic than the other types, since they can take into account such factors as temporal and spatial variations in meteorological fields and non-linear processes, such as complex chemical reactions. Unfortunately, due to problems associated with numerical ~eudo-diffusion (e.g. Long and Pepper, 1976) and limitations of speed and capacity of available computers, which result in insufficient spatial resolution, grid models are often unsatisfactory for urban- or regional-scale problems. Furthermore, in most urban models, sources are assumed to be located at some effective height computed independently, usually from one of the semi-empirical formulae summarized by Briggs (1970). These formulae arenot applicable when local wind circulations such as lake breezes become important. In summa~, Gaussian plume models cannot be used in unsteady and non-unifo~ flow conditions,
1933
c. M.
1934
SHEIH
grid models may be seriously handicapped by their limited spatial resolution, and the currently available plume-rise calculations are rather restrictive in their application. To cope with these problems, a puff model that includes the effects of wind shear and dynamic plume rise is proposed in the present paper.
2. THE MODEL
(a) General description The plume from a point, line or area source is approximated by a series of puffs emitted from the specific source. The location and the shape of each puff are determined by six tracer particles, represented by P1 for 1 = 1,2,..., 6 as shown in Fig. 1. These particles are released at the surface of the source and at the end points of orthogonal axes parallel to the coordinate axes, xt for i = l-3. The locations of the six particles are defined by eight spatial variables, C,; for i = l-8. The variables with two indices shown in the figure represent an average value, i.e. &, = (t, + e;,/z. The variables <, refer,to a fixed coordinate system and are used in keeping track of the absolute location of the tracer particles, while variables x1, x2 and xs, which refer to a coordinate system moving with the puff and are used in constructing the distribution of pollutant concentration around the pufT center. The eight variables <, are computed in each time step, allowing for advection by the mean wind, turbulent diffusion, wind shear and plume-rise entrainment. The concentration distribution for each pti is assumed to be Gaussian and is calculated by fitting an ellipsoid through the cluster of the six particles such that the half-lengths of the principal axes of the ellipsoid specify the values of the standard deviations (or dispersion coefficients). The concentration at a point of interest is obtained by summing thecontributions from nearby puffs. Further details of the model are given in the following sections. (b) Advection and da&ion For clarity in discussion, the movement of a puff due to advection and diffusion is explained in terms of a
Fig. 1. The coordinate system of tracer particles.
(d) 0
PLUME SEGMENT 51
19
w t
Fig. 2. Movement of tracer particles in a one-dimensional version of the puff model: (a) t = 0; (b) I = At; with advection only ; (c)I = At with advection and diffusion ; (d) at t = Twith r > At,(e) at T + At with advection only; and (f) at r + At with advection and diffusion.
one-dimensional model. Figure 2 shows the one dimensional spatial coordinate t, and the locations of the putt for various atmospheric conditions and at various times after release. The one-dimensional puff is represented by a thick solid line between the tracer particles P, and P2. At t=O, one tracer particle of the puff is placed at each end of the pollutant source. As the pollutant leaves the source, it will form a trace shown by the thick solid line. The amount of pollutant in each puff is the product of pollutant emission rate (q) and the time increment (At) used in segmenting the puff into a series of puffs. If the wind direction is toward the right and if there is no diffusion, the pollutant trace appears in the form of Fig. 2(b) at the end of the first time step (t = Art).The left boundary remains at the original location, whereas the right boundary has been advected for a distance of u2 At, where u1 is the wind velocity at P1. If there is diffusion, Fig. 2(c) shows that additional diffusive displacement A!, and At, will occur for the left and the right boundaries of the puff, respectively. It should be noted that, if the mean wind is strong enough, the puff cannot physically diffuse upwind on the first time step. This effect can be easily implemented in the present model but will not be considered here because the impact caused by such a simplification will not be significant on the simulation of the entire plume. After the first time step, each puff is detached from the source. Since the pollutant source no longer limits the movement of the puff, the displacement of the tracer particles becomes much simpler and consists of advection by the mean wind and expansion of the puff due to diffusion. The movement of the puff in this stage is illustrated in Fig. 2(d), (e) and (f). The above illustration can be easily extended to a threedimensional case and the locations of the six
A
particles at R + 1 time step can be written in terms of the parameters and variables at the previous time step n by the equations, for i = l-6,
r,n+k
=
t,+Si+(-
l)‘At,
& + uiAr -t (-
1931
puffpollutant dispersion model
for ?I =
0,
l)‘Afe,v for n 2 1,
(1)
between the two puff centers, i.e. C =2C’ at x, = u1 At, and x2 =x3 =O. Using Equations (2) and (3) under such condition, one can solve the half-width of the puff due to stretching of the puff by the mean wind during the first time step and obtain Q, = 0.4u, At.
where
This suggests that a better approximation variable SI is
Sf =
0,
{
u, At,
for upstream particles, for downstream particles,
u, v and ware the customary wind velocity com~n~ts in the directions of xi, x2 and x3; Pi represents the location af the ith tracer particle; and those variables without temporal superscript represent values at time step n. In short, the major three terms on the right hand side of Equation (1) are the requirement during the first time step that the upstream end of the puff has to remain attached to the source, the displacement of a tracer particle due to mean wind advection and turbulent diffusion. For simplicity in explanation, the previous illustration in Fig. 2 was assumed to be a top-hat or rectanguk distribution of pohutant concentration. In order to allow a series of discrete put% to approximate a continuous Gaussian plume, S, is determined by comparing the concentration distribution computed from the puff model with the corresponding Gaussian pIume. Let uI = a2 ; u1 = 0 for i = 3,. . ., 6; Af, = 0 for I= . 1, . . ., 6; 4 be emission rate of a continuous point source; and 0.5 A t be the time fraction travelled by the center of the puff, then the concentration C and C’ for the Gaussian plume and the puff can be written as C=
q 2x 62 CT3 u1
s, =
0.1 UsAt,
for upstream particles,
0.9 ut At,
for downstream particles.
(4)
The displacement due to turbulent diffusion can be estimated from the analytical solution givm by Sutton (1953). It suggests that the size of the puff along axis xf is e, = (2K,r)O5
(5)
where f, and K, are the puff width and eddy diffusivity, respectively, along xP It should be noted that many other formulae for the growth of plume width under various conditions are available (e.g. Gifford, 1977; and Sheih, 1977) and can be easily incorporated into the present model. Since the puff is followed in a Lagrangian frame, the turbulent properties can be considered to be relatively constant within a time step. A finitsdifkrence approximation of the time de rivative of Equation (5) is Aft =
K, OS5 2t At.
0
Since the horizontal wind velocity gradient is generally much smaller than the vertical, wind shear about the vertical axis will be neglected. Moreover, the wind shear will be assumed constant across each puff. With theseassumptions, only two additional variables, tt and 4s defined in Fig. I, are needed to determine the locations of tracer particles affected by wind shear. The two variables can be computed from the following expressions,
e:+”
= t, + 0.5 At [u(P,) - u{P,)f,
<‘B’&= s& + O.SAEL~[YIP~) - ~fPs)]_
where the coordinate origin is at the p&t source, and tag, b2 and b3 are the half-widths of the puff in the directions ofzcl, x2 and x3. A reasonable condition for estimation of the stretching of the puff by the mean wind during the first time step is to require pollutant concentrations computed from the plume equation and from summing the contributions from two adjacent puffs to be equal at the mid-point
for the
(7) (8)
The computation procedure to calculate the wind shear effect on the distribution of pollutant concentration is to find the principal axes of the sheared puff by fitting an ellipsoid through the six tracer parti& as discussed previously. It can be shown that the sheared puff can be expressed by the foilowing equation,
(9)
1936
c. hd. SHEIH
with u, = (tzt - &,_ &2,
for i = 1,2 and 3. (10)
The variables in Equation (10) are all single-indexed and should not be confused with the double-indexed variables defined previously. Thecorresponding distribution of pollutant concentration for the puff is C=
(IAt (2?c)3’2L,L,L, +(x2
exp(-
- $1)2.;i
OS[(x,
- zx1))0;2
(11)
+z+;$,
where L,, L, and L, are the principal axes of the ellipsoid in Equation (9). The formulae for computing these principal axes are L1 = (A cos2 a2 + Bsina,cosa,
+ Dsin2a2)-o*5
L2 = (ay2 sin2 aI + ag2 cos z al)-o.s t3 = (Asin’a,
Wf
011= arctan K&,); A = 0i2cOs2aI + ai2sina,; D = cq2[~:u;’
Fig. 4. Comparison of the concentrations computed from a series of sheared puffs and a Gaussian plume. At the plume center: 0 = Gaussian plume; and A = sheared puff model.
At one standard deviation from the plume center: 0 = Gaussian plume; v = above the ccntre; and 0 = below the center.
- Bsina2cosa2 + Dcos2a,)-0~5
with
B = 2~; 1 [g;‘<,
DOWNSTREAM DISTANCE(m)
cos a, + a;‘& sin a,] ; + &r;2
+ 11;
a2 = 0.5 arctan [B/(A - D)]. For clarity, a two-dimensional ellipsoid in Fig. 3 is used to illustrate the scales involved. The sheared puff model has been tested under the condition, u = 5 m s-l at the plume center, v = w = 0, h/ax, = aulax2 =o, dtdfax, =O.ls-l; Kr=5m2s-1 for i= 1,2, . . .. 6 and At = 10s. The results are compared with that of a Gaussian plume as shown in Fig. 4. At the center of the plume, the ~on~ntrations of the puff model and the Gaussian plume are almost exactly equal. However, in comparison with the Gaussian plume, the sheared puff model generates concentrations higher and lower at one standard de viation below and above the center of the plume, respectively. The departure from the Gaussian plume increases as the plume travels further downstream, due
to the accumulation of the shear effect. The situation can be visualized by imagining the plume cross sections perpendicular to the mean wind as a string of cards. Application of the previous shear rate is essentially pushing the upper and lower ends of the cards with faster and slower speeds, respectively. In comparison with the non-sheared case, the number ofcards per unit length, representing the concentration, will be larger below and smaller above, but there will be no change at the center line of the cards.
(d) P~@rise The governing equations for puff rise can be formulated by making assumptions similar to those of Morton et al. (1956). From the conservation laws of volume, rnorn~t~ and heat for an instantaneo~ source, the equations are i
;b3 (
= 4nb2j?jw,( + nb:1w,16(r),
(13)
>
$ ($nb3pw,) = $b3 ;
(0 - 0,)
+ &W&W,
(141
; CW3P(@ - W] = 4nb2B I WI,I P&d@ - @It) + ~a4.~o(ff. - @I?PO)* (13
Fig. 3. The coordinate system of the sheared putl.
where 6, 8, w, p. 6 and T are radius of the puff, entrainment constant, vertical velocity, density, potential temperature and temperature, respectively. The subscripts, o, a, R and p represent initial (exit) conditions, ambient conditions, reference conditions and the conditions of the puff relative to the ambient air. The variables without a subscript denote the con-
1937
A puff pollutant dispersion model
ditions inside the puff after leaving the source. The function s(t) is equal to unity for o c t < At and zero for all other values of t. Those terms with d(t) in the above equations represent the contribution (an addition to the original expression of Morton et al.) of the source to the pti before the model puff disengages from the stack. Since initial and ambient conditions are known, the unknown variables in the above equations are b, wP, 0 and p. However, if b can be predicted, the volume of the puff and the amount of ambient air entrained can be calculated and the density of the puff can be determined. Therefore, only three of the above four variables are independent. An iterative scheme is effective in solving the system of Equations (13)-(15) numerically because, if the radius of the puff is over-estimated, too much ambient air will be entrained and will result in reducing too much the temperature difference between the puff and the ambient air. Consequently, the buoyancy force and the vertical velocity will be reduced. This reduction in vertical velocity will reduce entrainment and thus will reduce the radius of the pti. In the derivation of the finite-difference equations for Equations (13)-( 15), the centered-time differences are used in derivatives and the average values of variables before and after an iteration are used for the non-derivative variables. The final problem to be resolved before computing the plume rise concerns the entrainment constant. Obviously, the entrainment constant developed by Morton et al. for a single puff or an instantaneous source is not applicable to the present model. This is due to the fact that in an urban model, sources are generally continuous and, with the present model, will generate a series of puffs overlapping each other. This will in turn reduce the amount of ambient air entrained into each of the puffs. To determine the entrainment constant for the present model, the puff-rise equations are solved for a continuous point source with various values of /? and fixed values for the rest of the parameters, namely b, = 1 m ; 0~ = TR= 8, = 290” C ; G = 0.05~~‘, At = 2 s; w, =Oms-‘; F, = 12.5m4sm2; and pn = 1.2 kgm-“, where the ambient vertical density gradient parameter, G = and the buoyancy parameter, F, - (gipR) +0x,, = 4nbz g (pa - p)/(3pR), are commonly used in plumerise studies (B&s, 1970). The results show that the best choice of the entrainment constant is /I = 0.14. Morton et al. (1956) found that the constant for an instantaneous puff is B = 0.285. This indicates that the effect of overlapping between successive puffs is essentially to block half of the surface area of each puff from entraining ambient air. In order to validate the present plume-rise model, an attempt is made to simulate the plume trajectory of Experiment No. 12 of the observations by Slawson and Csanady (1971). This case is chosen because the temperature lapse rate is near neutral and very uniform which enables us to assume a uniform eddy diffusivity for the simulation. The value for eddy diffusivity must be assumed, because information is
of 0
’
’ 1000
I
’
’
2ooo r/f,
’ 3ow
L
I
’
4ooo
Fig. 5. Comparison of modeled plume trajectories for eddy diffusivities of 0.5. 1.0 and 5.0m2s-’ (solid lines) and an observed plume trajectory (circles) obseki by S&son and Csanady (1971). The buoyancy length scale is e, = 0.43 m, and x and z are horizontal and vertical distances from the source, respectively.
insufficient for calculation. However, all other plumerise parameters such as stack-gas temperature, gas exit velocity, stack diameter, ambient wind velocity at the stack height and the atmospheric temperature profile are available in their report and are used in the simulation. Three values of atmospheric eddy diffusivity, 0.5, 1.0 and 5m2s-‘, are tested and the results are shown in Fig. 5. The buoyancy length scale (lb) used in normalizing the coordinates of the plume trajectory is 0.43m. An excellent agreement between the observed plume trajectory and the predicted trajectory occurs when the eddy diffusivity is 0.5m2 s-l. This value for a near-neutral stability condition appears to be reasonable (e.g. Tetsuji and Mellor, 1975). 3. PROCEDURES FOR GENERAL COMPUTATION
The techniques developed so far can be combined into a set of general computational procedures as follows : (1) variables at time step n are either initialized or predicted by use of their values at previous time step; (2) predict the locations of the tracer particles at time step n + 1, . . .. 6, t:;+I = r, +&+(-1)’ At
-, 2 for n = 0,
(16)
t;+'=[,+u,At+(-1)'
+B for n 2 1,
(l%I+Iw;+ll\)AI ^ \
L
/
(17)
1938
c. M. SWH
the grid system and treated with the finite-difference model.
(19) and estimate the equivalent puff radius b”+’ = o.s[(<;+’
Acknowledgements - The author would like to express his appreciation to Dr. W. J. Moroz, B. B. Hicks, G. D. Hess, and
J. D. Shannon for valuable discussions and reading of the manuscript.
- <;+I,
The time increment in Equation (16) is divided by 2 because the average exposure time of a puff in the first time step is only half of the time increment. (3) Solve Equations (13)-(15) for $+l and O”+‘, and iterate steps (2) and (3) until a stable value of b”+’ is obtained. (4) If concentrations are not needed at this time, update the variables, i.e. replacing variables at the nth time step by their values at the (n + 1)th step, and repeat the computations from steps (2)-(4). (5) If concentrations are needed, compute the principal axes with Equation (12) and concentrations with Equation (11) for each of the puffs and sum all the concentrations contributed by the nearby puffs to the particular point of interest.
REFERENCES Briggs G. A. (1970) Plume rise. AEC Critical Review Seria, USAEC Division of Technical Information Extension, Oak Ridge, Tennessee, U.S.A. Giflord F. A. (1977) Tropospheric relative diffusion observations. J. oppl. Met. 16, 311-313. Hilst G. R. (1967) An air pollution model of Connecticut. Proc. of IBM Scientific Computing Symposium on Water-Air Resource Management. Report No. IBM-320 1953, pp. 251-274. Lamb R. G. (1969) An air pollution model of Los Angeles. MS. Thesis, University of California at Los Angeles, U.S.A. Long P. E. and Pepper D. W. (1976) A comparison of six numerical schemes for calculating the advection of atmos-
pheric pollution. Preprint of the AMS Third Symposium on Atmospheric Turbulence, Diffusion and Air Quality, 19-22 Oc&r, Raleigh, N.C., pp. 181-187. Lucas D. H. 11958)The atmosnhexic nollution of cities. Int. J. Air Pol&
4. CONCLUSION
The unique feature of the present model is the incorporation of wind shear etTects and dynamic plume-rise calculation into a puff model. Application of the model in a typical atmospheric situation shows that wind shear is an important factor in determining pollutant concentration and should not be neglected. The results show that a more conventional model, which does not take into account the wind shear, overestimates and underestimates by a factor of 2 the pollutant concentrations at one standard deviation above and below the plume center, respectively, at 600 m downstream from the source. Comparison of a plume-rise simulation with an observation by Slawson and Csanady (1971) shows that excellent agreement on the plume trajectory can be achieved, if a proper atmospheric eddy diffusivity can be chosen. The present model should be particularly useful for modeling pollutant dispersion under highly unsteady and non-uniform flow conditions such as of heatisland, lake-breeze or sea-breeze circulations, where plume rise can not be accurately computed with conventional formulae. The present model is also recommended for use in treating s&grid-scale sources in a finite-difference grid model where grid dimensions are much larger than source dimensions. All soufca smaller than the grid increment could be treated initially with the present puffmodel and when the puffs have grown to a size comparable with the grid dimension, the pollutant in the puff could be passed to
1, 71-86.
1
a
Ludwig F. L., Gasiorek L. S. and Ruff R. E (1977) Simplification of a Gaussian puff model for real-time minicomputer use. Atmospheric Enoironmem 11.431-436. Morton B. R., Taylor G. I. and Turner J. S. (1956) Turbulent gravitational convection from maintained and instantaneous sources. Proc. R. Sot. A 234, l-23. Richards J. M. (1970) The effectof wind shear on a puff. Q. J. R. Met. Sot. %, 702-714. Roberts J. J., Croke E. S. and Kennedy A. S. (1970) An urban
atmospheric dispersion model. Proc. of Symposium on MuhidcSource Urban Diffusion Models. Air Poll. Cont. 08. GubIication No. AP-86, pp. 6.1-6.72. Schiermeier F. A. and Niemeyer L. E. (1970) Large power plant efIIuent study. National Air Poll. Cont. Administration Publication No. APTD 70-2, pp. 182-230. Sheih C. M. (1977) On the relative importance of singleparticle and relative diffusion for plume dispersion. I& mints of AMS Fiih Conference on Probabilitv and Statistics in Atmospheric Sciences, 15-18 Nov., Las &as, Nevada, pp. 265-268. Shieh L. J., Davidson B. and Friend J. P. (1970) A model of diffusion in urban atmosphere: SO1 in Greater New York. Proc. of Symposium on Multiple-Source Urban Diffusion Models. Air Poll. Cont. OK Publication No. AP-86, pp. 10.1-10.40. Slawson P. R. and Csanady G. T. (1971) The effect of atmosphaic conditions on plume rise. J. Fluid Me& 47, 33-49. Sutton 0. G. (1953) Micrometeorology, p. 134. McGraw-Hill, New York. Tetsuji Y. and Mellor (1975) A simulation of the Wangara atmospheric boundary layer data. J. atmos. Sci. 32, 2309-2329. Turner D. B. (1964)A diffusion model for urban area. J. appl. Met. 3, 83-91. Turner D. B. (1969) Workbook of Atmospheric Dispersion Estimates. U.S. Department of Health, Education and Welfare, National Air Pollution Control Administration, Cincinnati, Ohio, 84.
A PUFF POLLUTANT DISPERSION MODEL WITH WIND SHEAR AND DYNAMIC PLUME RISE* C. M. Atm~phe~c
SHEIH
Physics Section, Radiolo~~l and ~vironm~tal Research Division, Argonne National Laboratory, Argonne, IL 60439, U.S.A.
(First receiwd 8 Nuuember 1977 and in&al form 6 March 1978)
Abstract - A puff diffusion model, which includes wind shear and dynamic plume rise, is developed for numerical prediction of pollutant concentrations under unsteady and non-uniform flow conditions. The plume from a continuous source is treated as a series of put%emitted successively from the source. Each p&is represented by a set of six tracer particles, which define the size, shape and location of the puff. Initially these. particles are located at the surface of the source, on arbitrarily chosen orthogonal axes. The location of the particles is computed at each time step by taking into account advection, eddy diffusion, wind shear and entrainment of ambient air during plume rise. The concentration distribution of each puff isdetermined by fitting an ellipsoid to thecluster of the six particles and assuming a three-dimensional Gaussian distribution, with standard deviation quaf to the half-lengths of the principaf axes of the ellipsoid. The concentration at a point of interest is obtained by summing the ~n~butions from nearby puffs. The effect of wind shear on the pollutant con~tration is investigated by use of a typical wind shear encountered in the atmosphere. The resufts show that, at 6OOmdownstream from the source, the present mode1gives con~trations a factor of 2 higher and lower at one standard deviation below and above the plume center, respectiveiy, than that of conventional models in which no wind shear is considered. The phtme&e formulation is calibrated against the observations compiled by Briggs and the model is used to predict the trajectory of a plume observed by Slawson and Csanady. Excellent agreement between the prediction and the observation can be achieved if an appropriate eddy diffusivity is chosen.
1. iNTRODUCTION
A considerable number of diffusion models have been developed during the past several decades. Basically, the available models can be classified into plume, puiI and grid models. Plume models have been developed and applied by many investigators, such as Lucas (1958), Turner (19&Q, Hilst (1%7), and Shieh er al. (1970). In these models the growth of the plume width is a function of atmospheric stability and of travel time, as summarized by Turner (1969), and the concentration distribution of the plume in the cross-wind direction is assumed to be Gaussian. However, plume models are derived under the assumption of spatial and temporal uniformity in wind and atmospheric stability. Obviously, these assumptions are rather limiting for many air pollution problems in which factors such as non-~ifo~ to~graphical features, sea or lake breeze circulations, and urban heat-island effects can cause significant temporal and spatial variations in meteorological conditions. In order to include local effects, Lamb (1%9), Roberts er al. (1970) and Ludwig et al. (1977) have developed puff models. In these models, source emissions are treated as a series of puffs emitted into the atmosphere The concentration distribution within each pnfI is assumed to be Gaussian. In neither plume nor puff models is wind shear taken into account. This can be a rather serious * Work supported by U.S. Department of Energy and The Pennsylvania State University.
omission, as the observations by Schiermeier and Niemeyer (1970) show that plume cross-sections are highly elongated by wind shear. Although the e&t of wind shear on an isolated puff has been investigated by Richards (197(l), the results can not be applied directly to the present problem because the present puff is not an isolated one but a segment of a continuous plume. Moreover, Richards’ solution is obtained under constant shear and constant total buoyancy force, which is inappropriate for the non-uniform, unsteady state addressed in the present paper. In general, grid models solve finite-difference equations of diffusion. These models have the potential to be more realistic than the other types, since they can take into account such factors as temporal and spatial variations in meteorological fields and non-linear processes, such as complex chemical reactions. Unfortunately, due to problems associated with numerical ~eudo-diffusion (e.g. Long and Pepper, 1976) and limitations of speed and capacity of available computers, which result in insufficient spatial resolution, grid models are often unsatisfactory for urban- or regional-scale problems. Furthermore, in most urban models, sources are assumed to be located at some effective height computed independently, usually from one of the semi-empirical formulae summarized by Briggs (1970). These formulae arenot applicable when local wind circulations such as lake breezes become important. In summa~, Gaussian plume models cannot be used in unsteady and non-unifo~ flow conditions,
1933
c. M.
1934
SHEIH
grid models may be seriously handicapped by their limited spatial resolution, and the currently available plume-rise calculations are rather restrictive in their application. To cope with these problems, a puff model that includes the effects of wind shear and dynamic plume rise is proposed in the present paper.
2. THE MODEL
(a) General description The plume from a point, line or area source is approximated by a series of puffs emitted from the specific source. The location and the shape of each puff are determined by six tracer particles, represented by P1 for 1 = 1,2,..., 6 as shown in Fig. 1. These particles are released at the surface of the source and at the end points of orthogonal axes parallel to the coordinate axes, xt for i = l-3. The locations of the six particles are defined by eight spatial variables, C,; for i = l-8. The variables with two indices shown in the figure represent an average value, i.e. &, = (t, + e;,/z. The variables <, refer,to a fixed coordinate system and are used in keeping track of the absolute location of the tracer particles, while variables x1, x2 and xs, which refer to a coordinate system moving with the puff and are used in constructing the distribution of pollutant concentration around the pufT center. The eight variables <, are computed in each time step, allowing for advection by the mean wind, turbulent diffusion, wind shear and plume-rise entrainment. The concentration distribution for each pti is assumed to be Gaussian and is calculated by fitting an ellipsoid through the cluster of the six particles such that the half-lengths of the principal axes of the ellipsoid specify the values of the standard deviations (or dispersion coefficients). The concentration at a point of interest is obtained by summing thecontributions from nearby puffs. Further details of the model are given in the following sections. (b) Advection and da&ion For clarity in discussion, the movement of a puff due to advection and diffusion is explained in terms of a
Fig. 1. The coordinate system of tracer particles.
(d) 0
PLUME SEGMENT 51
19
w t
Fig. 2. Movement of tracer particles in a one-dimensional version of the puff model: (a) t = 0; (b) I = At; with advection only ; (c)I = At with advection and diffusion ; (d) at t = Twith r > At,(e) at T + At with advection only; and (f) at r + At with advection and diffusion.
one-dimensional model. Figure 2 shows the one dimensional spatial coordinate t, and the locations of the putt for various atmospheric conditions and at various times after release. The one-dimensional puff is represented by a thick solid line between the tracer particles P, and P2. At t=O, one tracer particle of the puff is placed at each end of the pollutant source. As the pollutant leaves the source, it will form a trace shown by the thick solid line. The amount of pollutant in each puff is the product of pollutant emission rate (q) and the time increment (At) used in segmenting the puff into a series of puffs. If the wind direction is toward the right and if there is no diffusion, the pollutant trace appears in the form of Fig. 2(b) at the end of the first time step (t = Art).The left boundary remains at the original location, whereas the right boundary has been advected for a distance of u2 At, where u1 is the wind velocity at P1. If there is diffusion, Fig. 2(c) shows that additional diffusive displacement A!, and At, will occur for the left and the right boundaries of the puff, respectively. It should be noted that, if the mean wind is strong enough, the puff cannot physically diffuse upwind on the first time step. This effect can be easily implemented in the present model but will not be considered here because the impact caused by such a simplification will not be significant on the simulation of the entire plume. After the first time step, each puff is detached from the source. Since the pollutant source no longer limits the movement of the puff, the displacement of the tracer particles becomes much simpler and consists of advection by the mean wind and expansion of the puff due to diffusion. The movement of the puff in this stage is illustrated in Fig. 2(d), (e) and (f). The above illustration can be easily extended to a threedimensional case and the locations of the six
A
particles at R + 1 time step can be written in terms of the parameters and variables at the previous time step n by the equations, for i = l-6,
r,n+k
=
t,+Si+(-
l)‘At,
& + uiAr -t (-
1931
puffpollutant dispersion model
for ?I =
0,
l)‘Afe,v for n 2 1,
(1)
between the two puff centers, i.e. C =2C’ at x, = u1 At, and x2 =x3 =O. Using Equations (2) and (3) under such condition, one can solve the half-width of the puff due to stretching of the puff by the mean wind during the first time step and obtain Q, = 0.4u, At.
where
This suggests that a better approximation variable SI is
Sf =
0,
{
u, At,
for upstream particles, for downstream particles,
u, v and ware the customary wind velocity com~n~ts in the directions of xi, x2 and x3; Pi represents the location af the ith tracer particle; and those variables without temporal superscript represent values at time step n. In short, the major three terms on the right hand side of Equation (1) are the requirement during the first time step that the upstream end of the puff has to remain attached to the source, the displacement of a tracer particle due to mean wind advection and turbulent diffusion. For simplicity in explanation, the previous illustration in Fig. 2 was assumed to be a top-hat or rectanguk distribution of pohutant concentration. In order to allow a series of discrete put% to approximate a continuous Gaussian plume, S, is determined by comparing the concentration distribution computed from the puff model with the corresponding Gaussian pIume. Let uI = a2 ; u1 = 0 for i = 3,. . ., 6; Af, = 0 for I= . 1, . . ., 6; 4 be emission rate of a continuous point source; and 0.5 A t be the time fraction travelled by the center of the puff, then the concentration C and C’ for the Gaussian plume and the puff can be written as C=
q 2x 62 CT3 u1
s, =
0.1 UsAt,
for upstream particles,
0.9 ut At,
for downstream particles.
(4)
The displacement due to turbulent diffusion can be estimated from the analytical solution givm by Sutton (1953). It suggests that the size of the puff along axis xf is e, = (2K,r)O5
(5)
where f, and K, are the puff width and eddy diffusivity, respectively, along xP It should be noted that many other formulae for the growth of plume width under various conditions are available (e.g. Gifford, 1977; and Sheih, 1977) and can be easily incorporated into the present model. Since the puff is followed in a Lagrangian frame, the turbulent properties can be considered to be relatively constant within a time step. A finitsdifkrence approximation of the time de rivative of Equation (5) is Aft =
K, OS5 2t At.
0
Since the horizontal wind velocity gradient is generally much smaller than the vertical, wind shear about the vertical axis will be neglected. Moreover, the wind shear will be assumed constant across each puff. With theseassumptions, only two additional variables, tt and 4s defined in Fig. I, are needed to determine the locations of tracer particles affected by wind shear. The two variables can be computed from the following expressions,
e:+”
= t, + 0.5 At [u(P,) - u{P,)f,
<‘B’&= s& + O.SAEL~[YIP~) - ~fPs)]_
where the coordinate origin is at the p&t source, and tag, b2 and b3 are the half-widths of the puff in the directions ofzcl, x2 and x3. A reasonable condition for estimation of the stretching of the puff by the mean wind during the first time step is to require pollutant concentrations computed from the plume equation and from summing the contributions from two adjacent puffs to be equal at the mid-point
for the
(7) (8)
The computation procedure to calculate the wind shear effect on the distribution of pollutant concentration is to find the principal axes of the sheared puff by fitting an ellipsoid through the six tracer parti& as discussed previously. It can be shown that the sheared puff can be expressed by the foilowing equation,
(9)
1936
c. hd. SHEIH
with u, = (tzt - &,_ &2,
for i = 1,2 and 3. (10)
The variables in Equation (10) are all single-indexed and should not be confused with the double-indexed variables defined previously. Thecorresponding distribution of pollutant concentration for the puff is C=
(IAt (2?c)3’2L,L,L, +(x2
exp(-
- $1)2.;i
OS[(x,
- zx1))0;2
(11)
+z+;$,
where L,, L, and L, are the principal axes of the ellipsoid in Equation (9). The formulae for computing these principal axes are L1 = (A cos2 a2 + Bsina,cosa,
+ Dsin2a2)-o*5
L2 = (ay2 sin2 aI + ag2 cos z al)-o.s t3 = (Asin’a,
Wf
011= arctan K&,); A = 0i2cOs2aI + ai2sina,; D = cq2[~:u;’
Fig. 4. Comparison of the concentrations computed from a series of sheared puffs and a Gaussian plume. At the plume center: 0 = Gaussian plume; and A = sheared puff model.
At one standard deviation from the plume center: 0 = Gaussian plume; v = above the ccntre; and 0 = below the center.
- Bsina2cosa2 + Dcos2a,)-0~5
with
B = 2~; 1 [g;‘<,
DOWNSTREAM DISTANCE(m)
cos a, + a;‘& sin a,] ; + &r;2
+ 11;
a2 = 0.5 arctan [B/(A - D)]. For clarity, a two-dimensional ellipsoid in Fig. 3 is used to illustrate the scales involved. The sheared puff model has been tested under the condition, u = 5 m s-l at the plume center, v = w = 0, h/ax, = aulax2 =o, dtdfax, =O.ls-l; Kr=5m2s-1 for i= 1,2, . . .. 6 and At = 10s. The results are compared with that of a Gaussian plume as shown in Fig. 4. At the center of the plume, the ~on~ntrations of the puff model and the Gaussian plume are almost exactly equal. However, in comparison with the Gaussian plume, the sheared puff model generates concentrations higher and lower at one standard de viation below and above the center of the plume, respectively. The departure from the Gaussian plume increases as the plume travels further downstream, due
to the accumulation of the shear effect. The situation can be visualized by imagining the plume cross sections perpendicular to the mean wind as a string of cards. Application of the previous shear rate is essentially pushing the upper and lower ends of the cards with faster and slower speeds, respectively. In comparison with the non-sheared case, the number ofcards per unit length, representing the concentration, will be larger below and smaller above, but there will be no change at the center line of the cards.
(d) P~@rise The governing equations for puff rise can be formulated by making assumptions similar to those of Morton et al. (1956). From the conservation laws of volume, rnorn~t~ and heat for an instantaneo~ source, the equations are i
;b3 (
= 4nb2j?jw,( + nb:1w,16(r),
(13)
>
$ ($nb3pw,) = $b3 ;
(0 - 0,)
+ &W&W,
(141
; CW3P(@ - W] = 4nb2B I WI,I P&d@ - @It) + ~a4.~o(ff. - @I?PO)* (13
Fig. 3. The coordinate system of the sheared putl.
where 6, 8, w, p. 6 and T are radius of the puff, entrainment constant, vertical velocity, density, potential temperature and temperature, respectively. The subscripts, o, a, R and p represent initial (exit) conditions, ambient conditions, reference conditions and the conditions of the puff relative to the ambient air. The variables without a subscript denote the con-
1937
A puff pollutant dispersion model
ditions inside the puff after leaving the source. The function s(t) is equal to unity for o c t < At and zero for all other values of t. Those terms with d(t) in the above equations represent the contribution (an addition to the original expression of Morton et al.) of the source to the pti before the model puff disengages from the stack. Since initial and ambient conditions are known, the unknown variables in the above equations are b, wP, 0 and p. However, if b can be predicted, the volume of the puff and the amount of ambient air entrained can be calculated and the density of the puff can be determined. Therefore, only three of the above four variables are independent. An iterative scheme is effective in solving the system of Equations (13)-(15) numerically because, if the radius of the puff is over-estimated, too much ambient air will be entrained and will result in reducing too much the temperature difference between the puff and the ambient air. Consequently, the buoyancy force and the vertical velocity will be reduced. This reduction in vertical velocity will reduce entrainment and thus will reduce the radius of the pti. In the derivation of the finite-difference equations for Equations (13)-( 15), the centered-time differences are used in derivatives and the average values of variables before and after an iteration are used for the non-derivative variables. The final problem to be resolved before computing the plume rise concerns the entrainment constant. Obviously, the entrainment constant developed by Morton et al. for a single puff or an instantaneous source is not applicable to the present model. This is due to the fact that in an urban model, sources are generally continuous and, with the present model, will generate a series of puffs overlapping each other. This will in turn reduce the amount of ambient air entrained into each of the puffs. To determine the entrainment constant for the present model, the puff-rise equations are solved for a continuous point source with various values of /? and fixed values for the rest of the parameters, namely b, = 1 m ; 0~ = TR= 8, = 290” C ; G = 0.05~~‘, At = 2 s; w, =Oms-‘; F, = 12.5m4sm2; and pn = 1.2 kgm-“, where the ambient vertical density gradient parameter, G = and the buoyancy parameter, F, - (gipR) +0x,, = 4nbz g (pa - p)/(3pR), are commonly used in plumerise studies (B&s, 1970). The results show that the best choice of the entrainment constant is /I = 0.14. Morton et al. (1956) found that the constant for an instantaneous puff is B = 0.285. This indicates that the effect of overlapping between successive puffs is essentially to block half of the surface area of each puff from entraining ambient air. In order to validate the present plume-rise model, an attempt is made to simulate the plume trajectory of Experiment No. 12 of the observations by Slawson and Csanady (1971). This case is chosen because the temperature lapse rate is near neutral and very uniform which enables us to assume a uniform eddy diffusivity for the simulation. The value for eddy diffusivity must be assumed, because information is
of 0
’
’ 1000
I
’
’
2ooo r/f,
’ 3ow
L
I
’
4ooo
Fig. 5. Comparison of modeled plume trajectories for eddy diffusivities of 0.5. 1.0 and 5.0m2s-’ (solid lines) and an observed plume trajectory (circles) obseki by S&son and Csanady (1971). The buoyancy length scale is e, = 0.43 m, and x and z are horizontal and vertical distances from the source, respectively.
insufficient for calculation. However, all other plumerise parameters such as stack-gas temperature, gas exit velocity, stack diameter, ambient wind velocity at the stack height and the atmospheric temperature profile are available in their report and are used in the simulation. Three values of atmospheric eddy diffusivity, 0.5, 1.0 and 5m2s-‘, are tested and the results are shown in Fig. 5. The buoyancy length scale (lb) used in normalizing the coordinates of the plume trajectory is 0.43m. An excellent agreement between the observed plume trajectory and the predicted trajectory occurs when the eddy diffusivity is 0.5m2 s-l. This value for a near-neutral stability condition appears to be reasonable (e.g. Tetsuji and Mellor, 1975). 3. PROCEDURES FOR GENERAL COMPUTATION
The techniques developed so far can be combined into a set of general computational procedures as follows : (1) variables at time step n are either initialized or predicted by use of their values at previous time step; (2) predict the locations of the tracer particles at time step n + 1, . . .. 6, t:;+I = r, +&+(-1)’ At
-, 2 for n = 0,
(16)
t;+'=[,+u,At+(-1)'
+B for n 2 1,
(l%I+Iw;+ll\)AI ^ \
L
/
(17)
1938
c. M. SWH
the grid system and treated with the finite-difference model.
(19) and estimate the equivalent puff radius b”+’ = o.s[(<;+’
Acknowledgements - The author would like to express his appreciation to Dr. W. J. Moroz, B. B. Hicks, G. D. Hess, and
J. D. Shannon for valuable discussions and reading of the manuscript.
- <;+I,
The time increment in Equation (16) is divided by 2 because the average exposure time of a puff in the first time step is only half of the time increment. (3) Solve Equations (13)-(15) for $+l and O”+‘, and iterate steps (2) and (3) until a stable value of b”+’ is obtained. (4) If concentrations are not needed at this time, update the variables, i.e. replacing variables at the nth time step by their values at the (n + 1)th step, and repeat the computations from steps (2)-(4). (5) If concentrations are needed, compute the principal axes with Equation (12) and concentrations with Equation (11) for each of the puffs and sum all the concentrations contributed by the nearby puffs to the particular point of interest.
REFERENCES Briggs G. A. (1970) Plume rise. AEC Critical Review Seria, USAEC Division of Technical Information Extension, Oak Ridge, Tennessee, U.S.A. Giflord F. A. (1977) Tropospheric relative diffusion observations. J. oppl. Met. 16, 311-313. Hilst G. R. (1967) An air pollution model of Connecticut. Proc. of IBM Scientific Computing Symposium on Water-Air Resource Management. Report No. IBM-320 1953, pp. 251-274. Lamb R. G. (1969) An air pollution model of Los Angeles. MS. Thesis, University of California at Los Angeles, U.S.A. Long P. E. and Pepper D. W. (1976) A comparison of six numerical schemes for calculating the advection of atmos-
pheric pollution. Preprint of the AMS Third Symposium on Atmospheric Turbulence, Diffusion and Air Quality, 19-22 Oc&r, Raleigh, N.C., pp. 181-187. Lucas D. H. 11958)The atmosnhexic nollution of cities. Int. J. Air Pol&
4. CONCLUSION
The unique feature of the present model is the incorporation of wind shear etTects and dynamic plume-rise calculation into a puff model. Application of the model in a typical atmospheric situation shows that wind shear is an important factor in determining pollutant concentration and should not be neglected. The results show that a more conventional model, which does not take into account the wind shear, overestimates and underestimates by a factor of 2 the pollutant concentrations at one standard deviation above and below the plume center, respectively, at 600 m downstream from the source. Comparison of a plume-rise simulation with an observation by Slawson and Csanady (1971) shows that excellent agreement on the plume trajectory can be achieved, if a proper atmospheric eddy diffusivity can be chosen. The present model should be particularly useful for modeling pollutant dispersion under highly unsteady and non-uniform flow conditions such as of heatisland, lake-breeze or sea-breeze circulations, where plume rise can not be accurately computed with conventional formulae. The present model is also recommended for use in treating s&grid-scale sources in a finite-difference grid model where grid dimensions are much larger than source dimensions. All soufca smaller than the grid increment could be treated initially with the present puffmodel and when the puffs have grown to a size comparable with the grid dimension, the pollutant in the puff could be passed to
1, 71-86.
1
a
Ludwig F. L., Gasiorek L. S. and Ruff R. E (1977) Simplification of a Gaussian puff model for real-time minicomputer use. Atmospheric Enoironmem 11.431-436. Morton B. R., Taylor G. I. and Turner J. S. (1956) Turbulent gravitational convection from maintained and instantaneous sources. Proc. R. Sot. A 234, l-23. Richards J. M. (1970) The effectof wind shear on a puff. Q. J. R. Met. Sot. %, 702-714. Roberts J. J., Croke E. S. and Kennedy A. S. (1970) An urban
atmospheric dispersion model. Proc. of Symposium on MuhidcSource Urban Diffusion Models. Air Poll. Cont. 08. GubIication No. AP-86, pp. 6.1-6.72. Schiermeier F. A. and Niemeyer L. E. (1970) Large power plant efIIuent study. National Air Poll. Cont. Administration Publication No. APTD 70-2, pp. 182-230. Sheih C. M. (1977) On the relative importance of singleparticle and relative diffusion for plume dispersion. I& mints of AMS Fiih Conference on Probabilitv and Statistics in Atmospheric Sciences, 15-18 Nov., Las &as, Nevada, pp. 265-268. Shieh L. J., Davidson B. and Friend J. P. (1970) A model of diffusion in urban atmosphere: SO1 in Greater New York. Proc. of Symposium on Multiple-Source Urban Diffusion Models. Air Poll. Cont. OK Publication No. AP-86, pp. 10.1-10.40. Slawson P. R. and Csanady G. T. (1971) The effect of atmosphaic conditions on plume rise. J. Fluid Me& 47, 33-49. Sutton 0. G. (1953) Micrometeorology, p. 134. McGraw-Hill, New York. Tetsuji Y. and Mellor (1975) A simulation of the Wangara atmospheric boundary layer data. J. atmos. Sci. 32, 2309-2329. Turner D. B. (1964)A diffusion model for urban area. J. appl. Met. 3, 83-91. Turner D. B. (1969) Workbook of Atmospheric Dispersion Estimates. U.S. Department of Health, Education and Welfare, National Air Pollution Control Administration, Cincinnati, Ohio, 84.