Environmental Fluid Mechanics, Lecture Notes- Physics, Study notes of Physics

Download Environmental Fluid Mechanics, Lecture Notes- Physics and more Study notes Physics in PDF only on Docsity! 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment 4A8: ENVIRONMENTAL FLUID MECHANICS Dispersion of Pollution in the Atmospheric Environment “Air and water, the two essential fluids on which all life depends, have become global garbage cans” “Mankind has probably done more damage to the Earth in the 20th century than in all of previous human history” (Jacques Yves Cousteau, 1910-1997) Prof. E. Mastorakos Hopkinson Lab Tel: 32690 E-mail: em257@eng.cam.ac.uk 1 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment Table of Contents 1. Introduction 2. Atmospheric structure, chemistry and pollution 3. Statistical description of turbulent mixing and reaction 4. Air Quality Modelling and Plume Dispersion 5. Turbulent reacting flows and stochastic simulations 6. Summary and main points 2 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment 2. Atmospheric structure, chemistry and pollution In this Chapter, we will give a quick revision of terms and concepts from chemical kinetics, which are needed to allow us to use and understand atmospheric chemistry. Some particular features of pollution chemistry then follow and information on the structure of the atmosphere is given. 2.1 Fundamental concepts – revision of chemistry Mole and mass fractions, concentrations Assume that pollutant A reacts with species B, which could be another pollutant or a background species (e.g. N2, O2, H2O in the air). The reaction rate depends on the amount of reactant present. There are many ways to quantify the amount of a species in a mixture: concentration, mole (or volume) and mass fractions are the most usual. The ratio of the number of kmols, ni, of a particular species i to the total number of kmols ntot in the mixture is the mole fraction or volume fraction: tot i i n n X  . (2.1) The mass fraction Yi is defined as the mass of i divided by the total mass. Using the obvious , (2.2) 1 11    N i i N i i YX where N is the total number of species in our mixture, the following can be easily derived for Yi and the mixture molecular weight MW : MW MW XY iii  , (2.3) 1 11             N i i i N i ii MW Y MWXMW . (2.4) Virtually always, the mixture molecular weight will be very close to that of air since the pollutant is dilute (i.e. even if it is heavy, its contribution to the weight of a kmol of mixture is very small). Equation (2.4) is included here only for completeness. The concentration (or molar concentration) of species i is defined as the number of kmols of the species per unit volume. The usual notation used for concentrations is Ci or the chemical symbol of the species in square brackets, e.g. [NO] for nitric oxide, or [A] for our generic pollutant A. From this definition and Eq. (2.1), V nX V n C totiii  , (2.5) and using the equation of state PV=ntotR 0T (R0 is the universal gas constant), we get: 5 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment TR P X PTRn nX C i tot toti i 00 /  (2.6) This relates the concentration to the mole fraction. In most atmospheric pollution problems, the concentrations are quoted in molecules/m3 or kmol/m3 and the volume fractions in parts per billion (ppb). A very common unit is kg of pollutant per m3 of air, which is the molar concentration times the molecular weight of the species. We can also relate the concentration to the mass fraction: i i i i i MW Y TR P MW MWY C   0 . (2.7) where  is the mixture density, e.g. the air density for our problems. Usually, the chemical reaction rate is expressed in terms of molar concentrations, while the conservation laws for mass and energy are expressed in terms of mass fractions. On the other hand, pollution monitoring equipment measures usually the volume fractions or the kg/m3 of the pollutant. The above relations are useful for performing transformations between the various quantities, which is very often needed in practice. Global and elementary reactions Chemical reactions occur when molecules of one species collide with molecules of another species and, for some of these collisions, one or more new molecules will be created. The chemical reaction essentially involves a re-distribution of how atoms are bonded together in the molecule. To achieve this, chemical bonds must be broken during the impact (i.e. the molecules must have sufficient kinetic energy) and other bonds must be formed. As the energy of these bonds depends on the nature of the atoms and on geometrical factors, the energy content of the products of the collision may be different from the energy content of the colliding molecules. This is the origin of the heat released (or absorbed) in chemical reactions. We write often that, for example, methane is oxidised according to CH4+2O2CO2+2H2O. This is an example of a global reaction. What we mean is that the overall process of oxidation uses 1 kmol of CH4 and 2 kmol of O2 to produce, if complete, 1 kmol of CO2 and 2 kmol of H2O. We do not mean that all this occurs during an actual molecular collision. This would be impossible to happen because it would involve too many bonds to break and too many bonds to form. However, the reactions CH4 + O  CH3 + OH O3 + NO  NO2 + O NO + OH + M  HNO2 + M are possible. For example, the first of these involves breaking one C-H bond and forming a O-H one. These reactions are examples of elementary reactions, i.e. reactions that can occur during a molecular collision. The overall chemical transformation follows hundreds or thousands of such elementary reactions and many species and radicals appear. By the term “radicals” we mean very reactive unstable molecules like O, H, OH, or CH3. The series of elementary reactions that describes the overall process is called a reaction mechanism or detailed chemical mechanism. Most chemical transformations occur following reaction mechanisms, rather than single reactions. The concept of global reaction helps us visualize the overall process and stoichiometry in an engineering sense. But when we identify the elementary reactions we can talk in detail about what really happens. 6 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment The Law of Mass Action A large part of the science of Chemical Kinetics is centred on identifying which elementary reactions are possible under various conditions for various species and to prescribe the rate, i.e. how quickly these reactions take place. This is given by the Law of Mass Action. Consider the generic elementary reaction MbPbPbMaRaRa MM  ...... 22112211 between reactants R1, R2, R3, …, from which products P1, P2, P3, …, are formed. M is an example of species that appears on both sides. The rates of reactants consumption and products formation and the reaction rate are given by: ... , ][ , ][ 2 2 1 1  a dt Rd a dt Rd  ... , ][ , ][ 2 2 1 1  b dt Pd b dt Pd   MM abdt Md  ][ (2.8) . (2.9) ...][][][ 321 321 aaa RRRk The parameter k is the reaction rate constant and [R1] is the concentration (in kmol/m 3) of reactant R1 etc. Equation (2.8) is a statement of the stoichiometry of the reaction: every a1 kmol of R1 is joined by a2 kmol of R2, etc., to produce simultaneously b1 kmol of P1, b2 kmol of P2 , etc.. If bM=aM, then M is called a third body: it may not be altered, but its presence is crucial for the success of the reaction, as it provides energy to, or takes energy away from, the collision between the reactants. Equation (2.9) is the Law of Mass Action and states that the reaction rate is proportional to the reactants concentrations, raised to their respective stoichiometric coefficients (i.e. a1, a2, etc). The amount of products does not affect . The reaction rate constant k is not a function of the reactants concentration and it is specific to the elementary reaction. Sometimes, Eq. (2.9) is used for a global reaction as an approximation. In that case, the reaction rate constant and the indices a1, a2, etc. are determined empirically. The reaction rate constant The reaction rate constant is given by the Arrhenius law:       TR E Ak act 0 exp (2.10) where A is the pre-exponential factor and Eact is the activation energy. These quantities come from experiment or statistical mechanics calculations. Out of all molecular collisions, only those with kinetic energy higher than the energy needed to break bonds inside the reactants’ molecules will result in reaction. The proportion of these collisions is given by exp(Eact/R0T) (from Kinetic Theory of Gases). The reaction rate constant increases very fast with temperature. In many environmental pollution problems, the flows are isothermal or the temperature changes little, so the reaction rate constant may be taken as uniform in space. However, the temperature may change significantly during the day so that large changes in reaction rates can be observed between noon and midnight. 7 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment VOC The term VOC refers to “volatile organic compounds”. By VOC we mean evaporated gasoline from petrol filling stations, vapours from refineries, organic solvents from paint and dry cleaning, and many others. These organic species are major contributors to smog and may be toxic or carcinogenic for humans. Their emissions are strictly regulated. Studies have found that glues and paints in our houses may have serious adverse health effects. There is a large effort today to switch to water-based paint, so that emissions of VOC’s are decreased. Heavy metals and dioxins If metals are contained in the reactants to our process (e.g. burner, kiln, chemical reactor), they may end up in the atmosphere (in pure form or in oxides) and they could be very dangerous (particularly Hg, Cd, Pb, As, Be, Cr, and Sb). This is a serious problem for municipal, toxic, and hospital waste incinerators. If the fuel contains chlorine, for example if it includes plastics (PVC), then there is a danger that dioxins may be formed  not at the flame, but on medium-temperature metal surfaces in the stack or inside the waste itself. Dioxins are chlorinated aromatic organic compounds whose chemistry is not very well known. They are extremely dangerous because they are carcinogenic even in concentrations of part per trillion. Municipal waste incinerators are needed to decrease the volume of waste going to landfills and to generate some power, but their use is controversial because of the danger of dioxins. If the plant operates at the design point, the exhaust is probably free of dangerous substances because a lot of attention has been given to the clean-up stage. However, open, uncontrolled burning of any waste or plastic material is extremely dangerous. Other sources of dioxins used to be the paper industry, where chlorine-containing compounds were used for bleaching the pulp and ended up in the waste stream of the paper mill, but this problem seems to have improved considerably with the introduction of new chlorine-free techniques. 2.3 Some chemistry of atmospheric pollution Smog formation One of the most important environmental problems is smog formation (Smoke + Fog), also called photochemical pollution. This is the brown colour we sometimes see above polluted cities and consists of a very large number of chemicals, but the most important are nitric oxide (NO), nitrogen dioxide (NO2), ozone (O3), and hydrocarbons (VOC). Some of the reactions participating in smog formation are discussed below. Oxidation of NO to NO2 Of the two possibilities: 2NO + O2  2NO2 R-I NO + O3  NO2 + O2 R-II the second is about four orders of magnitude faster. Therefore, the formation of nitrogen dioxide from the emitted NO depends on the availability of ozone. 10 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment Photolytic reactions A very important characteristic of smog is that it requires sunlight, which causes the photolysis of various chemicals into smaller fragments (‘photo’ = light; ‘lysis’ = breaking). Examples of sequences of reactions triggered or followed by photolysis: NO2 + hv  O + NO R-III O + O2 + M  O3 + M R-IV O + H2O  2OH R-V and NO + NO2 + H2O  2HNO2 R-VI HNO2 + hv  NO + OH R-VII The rate of photolytic reactions is usually given by an expression equivalent to Eq. (2.10), with the pre-exponential factor taken as a function of incident solar radiation. Hence the photolytic reaction rates are functions of the cloud cover and latitude in a given location, in addition to the time of the day and season. Smog is more common in sunny cities and in the summer because R-III (and other photolytic reactions) proceeds faster. The photostationary state Reactions R-II, R-III and R-IV suggest the following “cycle” for ozone: R-III produces an oxygen atom, which results in ozone formation through R-IV, which then reacts with NO in R-II. It is a common assumption that these reactions proceed very fast compared to others in smog chemistry and then the oxygen atom and ozone reach a quasi-steady state, called the “photostationary state”. If we assume that d[O]/dt =0 and d[O3]/dt=0 and we use the Law of Mass Action for the three reactions (R-II to R-IV), we get that ]NO[ ]NO[ ]O[ 23 II III k k  (2.11) where kII and kIII are the reaction rate constants for R-II and R-III respectively. Therefore we expect ozone (and smog) to be in high concentrations around midday when the photolytic reaction R-III has its peak. We also expect ozone to decrease to low values at night. Both these observations are approximately borne out by measurements. Above heavily polluted cities in the mornings, Eq. (2.11) is not very accurate due to the presence of VOC’s, which also produce NO2, as we shall see below. In the afternoons, when VOC’s tend to decrease due to their own photolysis, Eq. (2.11) is not a bad approximation. The hydrocarbons A prerequisite for smog is the presence of hydrocarbons, collectively called VOC’s. A global reaction describing their participation in smog formation is NO + VOC  NO2 + VOC* R-VIII where VOC* denotes some other organic radical. The NO2 will then follow the cycle in R-II, R-III and R-IV, which releases ozone. During the VOC chemistry, the major eye irritant CH3COO2NO2 may also appear (called peroxyacetyl nitrate or PAN). We see therefore that a combination of nitrogen oxide emission and hydrocarbons will under the action of sunlight cause serious pollution. The situation is worse above cities than above rural areas due to the high emissions from cars and 11 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment other activities. The chemistry of smog formation is extremely more complex than the above over- simplified picture and is the subject of intensive research. Acid rain Sulphur from oils and coal is usually oxidised during combustion and is then emitted in the form of SO2, which will react in the atmosphere according to: SO2 + OH  HSO3 R-VIII HSO3 + O2  SO3 + HO2 R-VIII SO3 + H2O  H2SO4 R-IX i.e. sulphuric acid is formed. This is absorbed on water droplets and may be deposited on the ground by rain. This is the notorious “acid rain”, a term that also includes nitric acid that is similarly formed. The situation has improved the last decade with the use of scrubbers in power stations, which capture SO2 before it is emitted. 2.4 Some comments on the atmosphere It is important to realize that the atmosphere at various heights behaves very differently, not only in terms of the chemistry that takes place, but also in terms of the fluid mechanics we see. Most pollution problems occur at the lower levels, although the ozone hole and global warming are issues of the higher levels too. Here, we present briefly the structure of the lower atmosphere, which will be necessary for understanding the atmospheric dispersion processes in Chapter 4. Atmospheric structure Pressure, temperature, “standard atmosphere” The pressure, density, and temperature of the atmosphere are related with height through g dz dp  (2.12) RTp  (2.13) pc g dz dT  (2.14) where R=R0/MWair. Equation (2.14) gives the Dry Adiabatic Lapse Rate (DALR) as 9.3 K/km. We very often use the “standard atmosphere” lapse rate, an average over all seasons of the year and across many regions of the globe, which is about 6.5 K/km. This can then be used in Eqs. (2.12) and (2.13) to find the vertical pressure and density distribution. The boundary layer and the troposphere, inversions The troposphere is the first 11 km above sea level and contains 75% of the mass of the whole atmosphere. It is usually divided into the boundary layer and the free troposphere. The boundary layer depth ranges between 500 and 3000 m and it is divided into the surface layer (the first 10% closest to the ground) and the neutral convective layer, also called the mixed layer. At the top of the boundary layer, we have the inversion layer. This structure is shown on Fig. 2.1. 12 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment Assume that the surface is uniform (i.e. same radiative properties and temperature everywhere) and let us neglect at this stage the presence of the atmosphere. Thermodynamic equilibrium implies that the net heat falling on the surface balances the heat emitted by the surface: R2(1-)S = 4R2Ts4 (2.15) giving Ts = [(1-)S / 4 ]1/4 (2.16) where Ts is the surface temperature,  the emissivity of the surface,  the albedo (the fraction of incident radiation that is reflected),  the Stefan-Boltzmann constant (equal to 5.67108 W/m2/K4). For Earth, the emissivity and the albedo are not constant (they depend on the nature of the surface, whether it is covered by water, ice, vegetation etc), but some average values are 0.61 and 0.3 respectively. These effective values include the presence of the atmosphere, clouds etc. Given that S  1367 W/m2, we get that Ts = 288 K. This is an estimate of the Earth’s average surface temperature, in the sense of an effective radiative temperature. Eq. (2.16) also shows the great sensitivity of the surface temperature in changes of solar radiation, albedo and emissivity; 1 K change can be brought about by 1-2% change in these parameters. Note that the greenhouse effect (discussed below) is already included in the effective values of emissivity and albedo. Area: R2 S Figure 2.3 Simplified view of Earth’s receiving radiation from the sun. In a more detailed view, Fig. 2.4 shows the various heat exchange processes taking place between the Earth’s surface and the incoming solar radiation (note that in this picture the solar radiation is expressed in terms of m2 of the total Earth’s surface area). The greenhouse gases (mainly water vapour, methane, and CO2) absorb some of the radiated heat from the surface and re-radiate it back to the surface, which results in a warming of the surface. If there were no greenhouse gases at all, 15 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment the surface would be very cold for human life, while if their concentration increases, the re-radiated fraction increases. The presence of atmospheric aerosols (natural or man-made) increases the reflection of the sun’s radiation, and therefore their presence can act so as to reduce the surface temperature. The power generation and the transport sectors produce significant amounts of aerosols and their effect on climate change is a very important topic of current research. Figure 2.4 Energy exchanges in the atmosphere (From the IPCC 4th Assessment Report; “Climate Change 2007: The Physical Basis”, CUP, 2007) 16 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment 2.6 Worked Examples Example 2.1 An air quality monitoring station measured one early morning that the volume fractions of NO, ozone and NO2 were 30, 40, and 170 ppb respectively. What were the corresponding concentrations in kmol/m3, in molecules/m3, and in kg/m3? The atmospheric conditions at the time of the measurement were 980 mbar pressure and 10ºC temperature. Solution From Eq. (2.6), , and using P=0.98x105 Pa, T=283 K, R0=8315 J/kmol/K, we get that 30x109 x 98x103 / (8315 x 283) = 1.24 x109 kmol/m3. Note that the 30 ppb becomes a volume fraction of 30x109. Similarly for the other species. To transform the kmol/m3 into molecules/m3, we need to multiply by Avogadro’s Number, 6.022x1026 molecules/kmol. To get the concentration in kg/m3, we multiply the molar concentration by the molecular weight of the species. TRPXC ii 0/ NOC Example 2.2 Estimate the mass of the air in the atmospheric boundary layer, assuming a total height of 3 km and a linear reduction of temperature with height of 6.5 K/km (the “standard atmosphere”) from a mean ground-level temperature of 288 K. Solution The standard atmosphere gives that zTzT 5.6)0()(  , where z is measured in km and T(0)=288K. The total mass of the air M(z) in the atmosphere up to a height z per unit area is given by . Writing and using the hydrostatic balance equation  z dzzM 0 )(  adzdT / gdzdP  / , differentiation of the equation of state RTp  gives: dzRdTdzRTddzdp ///   RagdzRTd   / azT dz R gaRd    )0(            )0( )0( ln )0( )( T azT aR aRgz       ln 1 )0( )0( )0()(          aR g T azT z  . So at z=3 km, the density of air is 0.9 kg/m3 and the temperature 268.5 K (4.5 C). The density variation can now be integrated to give:   1 , )()0( )0( )0( 1 1 )( 11 0     aR g mzTT aTm dzzM mm m z  . With a=6.5 K/km and z=3 km, the total mass per unit area becomes 3143.7 kg/m2. Note that if we had used the average value between [(0)+(3km)]/2 we would have obtained a mass per unit area of 3160.5 kg/m2. The fact that this estimate is so good reflects the fact that inside the boundary layer the density changes almost linearly with height (but not outside it). To find the total mass, we must multiply by the surface area of the Earth, which is 4Rearth2, and this gives a total mass in the boundary layer of 1.6x1018 kg. This corresponds to about 30% of the mass of the whole atmosphere. 17 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment Letting x go to zero, we obtain the species conservation equation: i ii w x m t Y       )( (3.2) Equation (3.2) is a partial differential equation (in time and space) and to be in a position to solve it, we need expressions for the mass flux and the rate of generation due to chemistry. The latter was covered in Section 2.1, while the former is discussed next. Mass flux, mass transfer and Fick’s Law of diffusion The mass flux im  for each species that appears in the species conservation equation is composed of two parts: an advective and a diffusive part. This result is given here without proof, as it can be proven from the Kinetic Theory of Gases (4A9, Part IIB). DIFFiADVii mmm ,,  (3.3) The advective mass flux is due to the bulk fluid motion and is given by: uYmYm iiADVi , (3.4) For the purposes of this course, the diffusive mass flux is given by Fick’s Law: x Y Dm iDIFFi    , (3.5) Fick’s Law states that the mass flux is proportional to the gradient of the mass fraction of the species. This is a diffusion process because it tends to make concentration gradients more uniform, i.e. it mixes the various species together. The coefficient D (m2/s) is the diffusion coefficient and, in general, depends on the nature of the diffusing species. For gases, it is a common approximation that the diffusion of heat and mass follow the same rate, i.e. D is related to the conductivity : pc D   (3.6) Equation (3.6) serves us to estimate D because tabulated values of conductivity and heat capacity are usually available. Throughout this course we assume that the diffusivity will be given by Eq. (3.6) with , , and cp taken as those of air at atmospheric conditions. Final instantaneous species conservation equation With these expressions, the species conservation equation takes the final form: i iii w x Y D xx uY t Y                   )()( (3.7) 20 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment It is important to know the physical mechanisms contributing to this equation: the first term in the l.h.s. corresponds to accumulation of species i, the second to advection by the bulk fluid motion, the first term in the r.h.s. corresponds to molecular diffusion and the last to the generation by the chemical reactions. In more dimensions and for a generic scalar  that is proportional to the mass fraction (e.g. our usual concentration in atmospheric pollution expressed in kg/m3), the governing transport conservation equation becomes: w x D x u t jj j          2 2 (3.8) in Cartesian tensor notation, where we have assumed an incompressible flow and a constant diffusivity, typically excellent assumptions in environmental fluid mechanics. If the scalar is inert, then simply . Equation (3.8) is our starting point for examining turbulent mixing in the following Sections. 0w 3.2 The averaged equations for a reactive scalar Averaged species conservation equation In a turbulent flow, we can write that the instantaneous mass fraction of a scalar is   and that the velocity is uuu  . It is easy to see that, by performing Reynolds decomposition and performing the averaging procedure (“The average of an average is the average”; “The average of a fluctuation is zero”; “The average of a product of fluctuations is not zero”) on Eq. (3.8), we get: w x D x u x u t jj j j j             2 2)()(  (3.9) The first term in the l.h.s. is the unsteady accumulation of , the second is due to mean advection, and the third is due to turbulent transport (or turbulent diffusion). The first term in the r.h.s. is due to molecular diffusion and the second is the mean reaction rate. Modelling the scalar flux – the eddy diffusivity It is usual engineering practice to model the turbulent transport term using the eddy diffusivity concept, also known as the gradient approximation. This model is motivated from the Kinetic Theory of Gases, where the mass flux is found to be proportional to the gradient of the mass fraction (Eq. 3.5) and the molecular diffusivity D is found to be proportional to the mean molecular speed and the mean free path between molecular collisions. By making an analogy between the random turbulent motions of “fluid particles” and the random molecular motion in a fluid, the turbulent transport term is written as j Tj x Du     (3.10) 21 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment with the eddy diffusivity DT given by (3.11) turbT LuCD  By a trial-and-error procedure and comparison with experimental data, the constant C is found to be around 0.1, but this depends on how Lturb is defined. There is a lot of criticism behind the use of the gradient approximation for modelling turbulent transport and indeed sometimes Eqs. (3.10) and/or (3.11) fail to predict the correct magnitude of  ju . Nevertheless, the eddy diffusivity concept remains a very useful approximation for providing a tractable closure to Eq. (3.9), which then becomes:   w x DD xx u t j T jj j                      (3.12) Note that DT may be a function of space and hence should be kept inside the derivative in the r.h.s. of Eq. (3.12). The eddy diffusivity concept is usually much better for an inert scalar than for a reacting scalar, but we use it anyway. For high Reynolds numbers TDD  , which suggests that the molecular diffusion may be neglected. To illustrate this, consider a wind flow of 5 m/s with a typical turbulence intensity of 10%, so that = 0.5 m/s. In the atmospheric boundary layer, the lengthscale is proportional to the height above the ground. Let us take that Lturb=500 m. Then DT = 25 m 2/s. At standard temperature and pressure, the molecular diffusivity of air is 2.2 x105 m2/s (convince yourselves with Eq. 3.6). Therefore the diffusion caused by molecular motions is negligible compared to the diffusion due to turbulence, which is a typical feature of turbulent flows at large Reynolds numbers. Molecular action is always present at the smallest (e.g. Kolmogorov) scales, but these contribute very little to the overall diffusion of the scalar (the small eddies just don’t “move far enough”). In other words, “where the smoke goes” is a function of the large scales only and the turbulent diffusivity suffices. u Governing equation for the fluctuations Starting from the instantaneous equation and performing the Reynolds decomposition, a series of transport equations for the higher moments may also be derived. Using   and uuu  in the instantaneous scalar equation (3.8) we get (before averaging):       ww x D x uu t jj jj          2 2)(  (3.13) Multiplying by   , expanding the derivatives and collecting terms gives:                           j j j j j j jj j x u x u x u t w x D x u t   2 2 w x D j     2 2 (3.14) 22 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment 3.3 The probability density function and averages Concept of pdf Figure 3.2 shows two possible time traces in a turbulent flow, e.g. of velocity u or of an inert or reactive scalar . If we measure how long the signal takes values between v and v+v, and then plot this quantity as a function of v, we get the curves on the right of the time trace. This is the probability density function of u. (v is called the random space variable of u; it is used, rather than u, for notational clarity.) In other words:  the probability of finding u in the region between v and v+v is P(v)v. t u(t) P(v) v (A) (B) Figure 3.2 Typical time series of signals with bimodal (A) or unimodal (B) pdfs. The dashed line shows the mean value. Comments on the pdf Cumulative distribution function: The probability of finding u<v is called the cumulative distribution function F(v). If the sample space of v is between a and b (i.e. the physics of the problem dictate that always a<u<b), then (3.20a)   v a tdtPvF )()( 0 (3.20b) )( aF 1 (3.20c) )( bF The probability of finding u between v1 and v2 is 25 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment (3.20d)   2 1 )()()( 12 v v tdtPvFvF which implies that dv vdF vP )( )(  (3.20e) Normalization condition: (3.21)   b a dvvP 1)( Mean and variance For any function y = f(u), if u is distributed according to P(u) (for simplicity, we now drop the difference between the sample space variable v and the real random variable u), we have that the average of y is given by:  b a duuPufufy )()()( (3.22) The average value of u itself is  b a duuuPu )( (3.23) In turbulence, it is convenient to work in terms of the fluctuation uuu  . Then, the mean of the fluctuation is 0)()()()(   uuduuPuduuuPduuPuuuuu b a b a b a (3.24) So the mean of the fluctuation does not tell us anything useful. We would particularly like to know how large is the excursion around the mean (i.e. how “far” the signal in Fig. 3.2 travels away from the mean). For this we use the variance, which is defined as   b a b a duuPuuduuPuu )()()( 2222 (3.25) The variance gives us a measure of how “wide” the pdf is. The quantity  (i.e. the square root of the variance) is called the root mean square (rms) of u. So, the pdf in Fig. 3.2A has a higher rms than the pdf of Fig. 3.2B, even if the mean is the same. It is clear that the mean alone does not convey all the information about the signal: the behaviour of Fig. 3.2A is very different from the behaviour of Fig. 3.2B. 26 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment The source term problem revisited (a) Let  kwy  )( , i.e. a first-order reaction. Then: .)()()()(  kdPkdPwwy b a b a    (b) Let , i.e. a second-order reaction. Then: 2)(  kwy     b a b a b a dPkdPkdPwwy  )()()()()()( 22 )()()2( 222222    kkkdPk b a Therefore the mean reaction rate of a second-order reaction depends not only on the mean, but on the variance of the scalar as well (equivalently, on the width of the pdf). We recovered Eq. (3.19) for A=B. Important conclusion for averages If a function is non-linear, the mean of the function is not equal to the function of the mean. For the special case of the linear function, the mean of the function is always equal to the function of the mean. In other words, if , then )(xfy  )()( xfxfy  . The importance of the fluctuations Assume that the trace of Fig. 3.2 corresponds to the concentration of the pollutant SO2 at a house located downwind of a chemical factory with aged scrubbers. Assume that the “safe” level of inhaling the pollutant is exactly at the mean value. Would you rather breathe from a plume resulting in a concentration following the curve of Fig. 3.2A or of Fig. 3.2B? Clearly, when the fluctuations are large, even if the mean value is deemed “safe”, the receptor is exposed often to high dosages of the pollutant. If the health effects after exposure are a very non-linear function of the pollutant concentration, then a regulation expressed only in terms of the mean value of the pollutant is not safe enough. This highlights the importance of Eq. (3.17), which gives the variance of the pollutant fluctuations. Ideally, the whole pdf of pollution concentration should be considered in conjunction to the “threshold” level for safety. However, this is not done often, not least because there are no easy tools available to predict the whole pdf for a reacting scalar. A way that this can be achieved is through a Monte Carlo (stochastic) simulation and this will be presented in Chapter 5. 3.4 Multi-variate probability density functions Fundamental properties The concepts above are readily generalised to many variables. Let 1 and 2 be random variables (e.g. turbulent reacting scalars) with random space variables 1 and 2, (for the sake of 27 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment Jointly-Gaussian pdf A usual shape for joint pdf’s in turbulence is the jointly Gaussian or joint normal pdf (Fig. 3.4):    2 1221 21 12 1 ),(  P                        2 2 2 22 21 221112 2 1 2 11 2 12 2 )())(( 2 )( )1( 1 exp        (3.31) We will encounter this equation again in Chapter 4 in the case of plume diffusion. v v 1 v 1 2 v 2 2 1  12 v v v 1 v 2 2 1  12  12 Figure 3.4 The joint normal probability density function for uncorrelated, positively and negatively correlated variables. Equal-probability contours are shown. Uniform distribution If  is equally probable to take any value between a and b, then P() is called the uniform distribution and the pdf is simply equal to )/(1)( abP  . In two dimensions and for uncorrelated variables, )( 1 )( 1 ),( 2211 21 abab P   . This situation is not very common in turbulence for a 30 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment scalar or the velocity. However, the uniform distribution is extremely useful as an analytical tool because other distributions can be related to it. Bounded pdfs In many situations of interest, the scalar is bounded. For example, if  corresponds to a concentration of pollutant in air, it cannot be less than zero nor higher than the maximum value it had at the source. In this case, a=0 and b=max. The pdf should, strictly speaking, be represented by an expression that reflects this bounded character. Various techniques to achieve this have been proposed. One is to use a “clipped Gaussian distribution”, where delta functions of variable strengths are used at the extrema points (a and b). Another is to use the so-called “beta function pdf”. Both have been used in turbulent reacting flows, but are analytically complex to represent. Qualitative shapes of pdfs for inert scalars corresponding to a typical mixing flow are given in Fig. 3.5. For a deeper discussion, see Bilger (1980). Pdfs of reactive scalars Unfortunately, the chemical reaction changes the shape of the pdf in an unknown way. Hence, no simple shape can be provided. The situation is even worse for the joint pdf of two reactive scalars. However, there are methods by which the joint scalar pdf can be numerically evaluated. It is evident from Eqs. (3.19) or (3.28) that knowledge of the shape of the joint pdf is sufficient to close the turbulent reacting flow problem, and hence a lot of research emphasis has been put in developing such methods (Pope, 2000). We will present a very simplified but still powerful technique for estimating the statistics of random phenomena in Chapter 5. 31 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment    A B 0 1   A 1 1 2 2 33 1 2 3 0 1   B1 2 3 Figure 3.5 Pdf shapes of a bounded inert scalar for a typical turbulent mixing flow. 32 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment H W L x y z U Figure 4.1 A typical city box model. By convention, we take the x-direction aligned with the wind. Further applications Emphasis on chemistry Equation (4.2) is not restricted to a single pollutant. Various researchers use it for examining also the chemistry and emphasize the reaction rate term, rather than the wind transport. In such box models, the governing equation for each species i becomes i i w dt dc  (4.4) Equation (4.4) is solved by numerical methods subject to a particular set of initial conditions and the solutions can help identify how the various pollutants are transformed during the day. Comparisons with experimental data can then assist in developing chemical mechanisms, such as “tuning” the rates of the various reactions like the smog-forming reactions R-II to R-IV in Chapter 2. Developing detailed chemical mechanisms for atmospheric pollution is a very active research area at present, not least because more and more chemicals come under regulation. Emphasis on inhomogeneity It is not a bad assumption to take the pollutant concentration as uniform in the z-direction, especially during unstable stratification. However, the assumption of homogeneous concentration in 35 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment the wind direction is usually much worse because often q is a function of x. This can be partly dealt with by re-deriving Eq. (4.2) for a thin strip of thickness x and hence obtaining a differential equation for dc/dx. (The derivation is given as an exercise in the Examples Paper. See also Example 4.2.) Emphasis on yearly averages Very common in Air Quality Modelling, Eq. (4.3) is used for a range of wind directions and speeds and a range of mixing heights, so that various meteorological conditions can be examined to find the corresponding pollutant concentrations. These are then weighted by the probability of occurrence of these particular conditions and hence a yearly average pollutant concentration can be calculated. Such calculations are important, e.g. for planning anti-pollution measures, for calculating the extra environmental burden of new industrial plants, etc. 4.2 Gaussian dispersion models Model problem The paradigm problem concerning pollution relatively close (e.g. a few km) to a source is the “chimney plume”. This is shown on Fig. 4.2. We are interested in: (a) the width of the plume downwind; (b) the concentration of the pollutant across the plume and particularly on the ground; (c) the difference between a steady emission (a “plume”) and an unsteady emission (a “puff”). To calculate these items is very important from a practical point of view and forms the topic of this Section. The material in this Section is based mostly on De Nevers (1995) and Csanady (1973). x y z U 0 z Figure 4.2 A typical chimney plume, showing the definition of the quantities used in the analysis. 36 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment Turbulent diffusion of an unsteady puff Derivation In many cases, we are interested in the way a pollutant spreads under the action of turbulent diffusion, when the emission of the pollutant is not continuous, but occurs only for a short time. Our full governing equation (Eq. 3.12) is our starting point, but to demonstrate how the solution comes about, let us assume zero mean velocities and finite and spatially-uniform turbulent diffusivity. Then, Eq. (3.12) becomes: w x K t j j       2 2 (4.5) The overbar denotes the mean concentration of the pollutant and the eddy diffusivity is now allowed to be a function of the direction. This is usually the case in atmospheric turbulence: the velocity fluctuations are not really equal in the three directions, and hence the eddy viscosity is strictly speaking a (diagonal) tensor. Hence the appearance of the subscript j in K in Eq. (4.5). We have also denoted the turbulent diffusivity by the symbol K (rather than DT) to conform to the standard notation in atmospheric pollution. For an inert pollutant, Eq. (4.5) is identical to the unsteady heat conduction equation for an “instantaneous source” and has a known solution. If Q kg/s of pollutant are released over a (very short) time t at point (x0,y0,z0), the solution of Eq. (4.5) gives for the mean pollutant concentration  in one, two, and three dimensions (in kg/m, kg/m2, kg/m3 respectively):             xx K xx tKt tQ x 2 0 2/12/1 )( 4 1 exp )()(2 )(   (4.6)                       yxyx K yy K xx tKKt tQ yx 2 0 2 0 2/1 )()( 4 1 exp ))((4 ),(   (4.7)                         zyxzyx K zz K yy K xx tKKKt tQ zyx 2 0 2 0 2 0 2/12/3 )()()( 4 1 exp )()(8 ),,(   (4.8) where t is the time from the release. Note the “symmetry” of the terms in the exponential, but also note that the behaviour of the maximum concentration at the centre of the cloud (i.e. at x0,y0,z0) has a different scaling with time depending on the dimensionality of the problem. In practical atmospheric dispersion of pollutant clouds, the wind has to be taken into account. The situation is visualized in Fig. 4.3. The governing equation becomes (for an inert pollutant) 2 2 2 2 2 2 z K y K x K x U t zyx                (4.9) The solution of Eq. (4.9) is again Eq. (4.8), but now we must interpret the time t as the downwind distance X of the centre of the pollutant cloud, divided by the wind speed U, and x,y,z as the distances from the centre of the cloud. In Eulerian coordinates (more useful!), x should be replaced 37 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment where u , , are the r.m.s. turbulent velocities in the three directions, assumed constant. Equation (4.13) can then be integrated (after a very long and difficult procedure) to give:  v w                                      ru Ux erfc ru xU ru Ux u U rwv uQ zyx 2 1 2 exp 2 1 2 exp )2( ),,( 22 22 2 2 22/3    (4.15) with 202 2 2 2 2 22 )( zz w u y v u xr        . Usually, the plume is slender, which implies that r  x , and usually . Both these assumptions are used to simplify Eq. (4.15) into: 1)1.0(/  oUu                     2 2 0 2 2 22 exp 1 2 ),,( zyzy zzy U Q zyx   (4.16) This is our final result for the plume diffusion problem and is one of the most important equations in this set of Lectures. It is often called the “Gaussian plume equation”. The x- dependence comes indirectly through the dependence of the disperion coefficients on x, while the dependence on y and z comes directly from the exponential. Diffusion from a point source The slender plume approximation is equivalent to neglecting the axial diffusion because, compared to the advection, it causes little spreading. Equation (4.16) is also the solution of the general governing equation for a statistically-steady inert scalar being injected at a rate Q at point (0,0,z0), with x-wise diffusion neglected: 2 2 2 2 z K y K x U zy          This problem is also known as “diffusion from a point source” and its solution is Eq. (4.16), with dx dU K yy 2 2   and dx dU K zz 2 2   . Experimentally, the K’s are found indirectly by measuring the concentration across the plume, applying Eq. (4.16) to the measured profiles, and hence evaluating  . This has proven to be a very useful tool to understand atmospheric turbulence and dispersion. The eddy diffusivities so found are called “apparent eddy diffusivities”. Short and long time diffusion For short distances, we have that , UxvtvK y / 22  UxwtwKz / 22  zzLz LwTw  , 2 , while for long distances from the source we obtain , , where TL,y is the Lagrangian autocorrelation time in the y-direction, Ly is an equivalent lengthscale, and similarly for the z-direction. Because in general yyLy LvTv  , 2 K w K v  and zLyL TT ,,  , in general zy KK  . Equation (4.12) for the cloud and Eq. (4.16) for the plume turn out to be valid for both short and long time mixing, with the difference that the dispersion coefficients are a different function of x in the two cases. From Eqs. (4.11) and (4.14), the width of the plume increases linearly with x for short times, while it increases as x1/2 for long times. 40 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment The difference between short and long-time diffusion is important for the theoretical study of turbulent mixing. If the release occurs at a length of the order of the Kolmogorov scale, then there is an additional, “very early” time regime, where the plume grows by the action of molecular motions, and there too, it grows as x1/2, but with a molecular, rather than a turbulent, diffusivity appearing in the proportionality constant. This is irrelevant to atmospheric pollutant problems, as the emission lengthscale – the stack diameter – is much larger than the Kolmogorov scale. The Gaussian Equation (4.16) is identical to a joint Gaussian expression in Section 3.5, Eq. (3.31), with 12=0. This is no coincidence: in order to derive Eq. (4.16) we have implicitly assumed that the turbulent motions in the two directions are uncorrelated. This would not be so if there were shear in the air flow, which would give rise to finite Reynolds stresses. The presence of Reynolds stresses would make a cross-section of the plume look skewed (c.f. Fig. 3.4) and Eq. (4.16) would need correction. This is often done for predicting dispersion in the atmospheric boundary layer, but is beyond the purposes of this course. Self-similar solutions There is an alternative way to derive Eq. (4.16). This is based on the concept of self- similarity, which is presented (for one-dimensional spreading) in the Appendix. Measurements of the dispersion coefficients In order to use Eq. (4.16) in practical applications, we need to know the dispersion coefficients y and z. These will be functions of downwind distance x and will also be unequal, since the characteristic turbulence velocities are expected to be different in the horizontal and vertical directions. Finally, we expect them to depend on atmospheric conditions: for stable conditions, the r.m.s. velocities are very low and hence the dispersion coefficients small; for unstable conditions, the turbulence is vigorous and we expect much larger dispersion coefficients. Measurements over a range of stability conditions give the results plotted in Fig. 4.4. Due to the complexity of atmospheric turbulence, the dispersion coefficients increase with x at a rate between 1 and 1/2, with the vertical dispersion coefficient increasing at a variable rate with x and the horizontal at about 0.9. The various curves correspond to various stability classifications, described in Table 4.1. From this Table and Fig. 4.4, the ’s can be estimated at a given downwind distance and hence Eq. (4.16) can be applied to give the mean pollutant. Table 4.1 Stability categories of Pasquill (adapted from De Nevers, 1995). A to F refer progressively from the very unstable to the very stable conditions, with D the neutral. Day Solar intensity Night Wind speed at 10 m high (m/s) Strong Moderate Slight Overcast Clear 0-2 A A-B B - - 2-3 A-B B C E F 3-5 B B-C C D E 5-6 C C-D D D D > 6 C D D D D 41 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment (a) Figure 4.4 Horizontal (a) and vertical (b) dispersion coefficients as a function of downwind distance from the source for the “standard” stability conditions of Table 4.1 (De Nevers, 1995, quoting Turner). 42 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment with the x-dependence coming from the x-dependence of the dispersion coefficient (e.g. Fig. 4.4). Note that we still need an image line source, but now both sources coincide at z0=0 (hence the 2 in the numerator). Equation (4.20) for the plume with vertical homogeneity is equivalent to diffusion from a vertical line source of strength Q/H, of course without any reflection since no boundaries exist normal to the y-direction. Plume rise Very often, the exhaust gases from chimneys are hotter than the environment or they have substantial momentum. This makes them rise and hence the release height z0 that should be used in the Gaussian puff and plume equations does not correspond truly to the stack height hs. Instead, it corresponds to the stack height plus the “plume rise” distance hpr, which is the approximate height above the stack before the plume turns 90 to follow the wind. To predict buoyant plume rise theoretically would involve a combined examination of momentum and buoyancy in the curved plume and is a little involved. An empirical relationship used sometimes by engineers designing stacks is:         s atmss pr T TT PD U DV h 68.25.1 (4.22) with hpr the plume rise (m); Vs the stack exit velocity (m/s); D the stack diameter (m); P the pressure of the atmosphere (bar); and Ts and Tatm the temperatures of the stack gas and the atmosphere respectively (K). 4.3 Effects of local meteorology Effect of inversions We mentioned in Section 4.2 that the case of one-dimensional diffusion in the horizontal direction approximately corresponds to the case where the plume is perfectly mixed in the vertical direction. If the mixing height is higher than the source height, i.e. if H > (hs+hpr), then there is every possibility that the plume will hit the inversion lid from below, from where it will be reflected downward etc. This situation is shown in Fig. 4.5, which also demonstrates how the image source idea can be used to account for the mixing height. It is clear that the vertical direction will eventually become quite homogeneous. 45 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment x z U z 0 H Image source Real source Figure 4.5 Effect of inversion lids on diffusion from a source below the inversion. Equation (4.20) applies in this case for distances far enough from the source. Since vertical diffusion is eliminated, the ground concentration decays as x1 (since, approximately, y ~ x1; Fig. 4.4), in contrast to a faster decay if the diffusion were two-dimensional (in the 2-D case, from Eq. 4.19 the maximum ground concentration decays by an additional z1). Therefore, low mixing heights will result in higher ground concentrations. Inversions assist in “trapping” pollution above a city, which can be thought as composed of a very large number of point sources. Cities with serious smog problems like Los Angeles, Athens, Istanbul and others, are located in sunny regions, bounded by mountains, and are close to the sea. In addition to more intense photochemistry, these conditions are favourable for the creation of inversion lids and hence such cities encounter high levels of pollution often. Plume shapes Depending on the conditions close to the plume source, the plume may acquire shapes very different from the “regular” diffusion shape defined by Eqs. (4.16) (without ground effect) and (4.17) (with ground effect). Various plume shapes and associated phenomena are described next and shown in Fig. 4.6. Fumigation If a plume is emitted into a very stable atmosphere (e.g. during a clear night), the plume travels downwind with very little dispersion. After dawn, as the sun heats the ground, an unstable layer is formed and the mixing height associated with the edge of this layer travels upwards as time 46 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment progresses. When this turbulent layer hits the plume, it causes very quick transport of the pollutant to the ground. This episode is called fumigation and it can result in very high ground concentrations, but for relatively short times. Lofting If either due to a very high stack or due to very high buoyant plume rise the plume source is effectively above the inversion, then the plume is dispersing upwards, but not downwards. The inversion acts like a solid boundary. This situation protects the ground receptor from high pollutant levels. Fanning If the plume source is within an inversion layer, it spreads horizontally but not vertically. This is called fanning and in such situations, the plume can travel long distances with little dilution. Trapped plume If the plume is below the inversion, it cannot penetrate it and the pollutant is not dispersing upwards. This causes higher ground concentrations than if the inversion were not there and is often the reason for high pollution levels in cities. Downwash Often it may happen that the stack sheds vortices or that the stack is in the wake of a large building. In both cases large eddies may trap the plume and push it towards the ground. These situations are called downwash or flagging and downdraught respectively. A rule of thumb by engineers is to use a stack height at least 2.5 times the height of the nearest building. Good aerodynamic design of the stack can help in reducing downwash. The same effect appears when the source is situated near a mountain, where the katabatic winds may not allow normal, Gaussian dispersion, but cause a bulk movement of the plume towards the ground. Katabatic winds also favour the creation of inversion lids by bringing cooler air close to the ground. Looping In very unstable situations, the eddies causing plume dispersion may be much larger than the plume width. This results in a “solid-body” motion of the plume, also called the “meandering plume”. See Section 5.3 for more discussion on this. Thermalling If the plume is very buoyant, it may break up in distinct parcels, as the plume itself is initiating thermals. 47 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment 4.5 Worked examples Example 4.1 A box model for a city has W=5 km, L=15 km, the mixing height is H=1 km, and the wind is blowing along the long dimension of the city with U=5 m/s. The background concentration of CO coming into the city is b=5x109 kg/m3, and the emission rate is q=4x109 kg/sm2. Find the concentration of CO above the city. Solution This is a straightforward application of Eq. (4.3) and is intended to give an idea on the order of magnitude of the quantities used in conjunction with city box models. Putting numbers gives cbm=25x10 9 kg/m3. Note that the result is independent of the width of the city. Example 4.2 The worst assumption in the standard box model is that the concentration is uniformly distributed along the wind direction. By doing a mass balance on infinitesimally thin strips normal to the wind, the following o.d.e. can be derived for c(x) (the derivation is given as an exercise in the Examples Paper): . Assume now that the pollutant undergoes photolysis, so that it is being destructed at a rate wHxqdxUdc  /)(/ kcw  . Find c(x) if q is constant and if no pollutant enters the city. Solution Since no pollutant enters the city, c(0)=0. Then, solving the o.d.e. gives: dx c U k UH q dc c U k UH q dx dc w H xq dx dc U     )(      )exp(1)( x U k kH q xc Therefore, the concentration of c will become equal to q/kH at x long compared to U/k. If photolysis were absent (k=0), integrating the o.d.e. would give that c=qx/UH, i.e. the pollutant is continuously increasing with downwind distance and we recover Eq. (4.3) for x=L. Photolysis reactions and other removal mechanisms (e.g. deposition to ground by rain) are important in controlling the levels of atmospheric pollution. Example 4.3 An accident in a chemical factory releases 1000 kg of SO2 from a chimney. The pollutant cloud is assumed to rise vertically upwards very quickly and then it is transported horizontally by the wind blowing with 5 m/s. The three turbulent diffusivities are estimated as 100 m2/s. Estimate the width of the cloud 5 km downstream and the pollutant concentration at the centre of the cloud. Solution This is a straightforward application of Eq. (4.10). The concentration at the centre will be: 50 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment 2/12/30 )()/(8 ),0,( zyx KKKUX q zX    , with q=1000 kg, X=5 km, Kx=Ky=Kz=100 m2/s, and U=5 m/s. The width of the cloud is estimated as UKX /4 (since we want the total width). Putting numbers, we get that 2SO 7.1x10 7 kg/m3 and the cloud width is 1265 m. Example 4.4 The effluent stream in a refinery is at 1 bar and 300 K and flows at 10 m/s from a chimney of diameter of 200 mm. The exhaust gases contain 100 ppm (by mass) of VOC’s and the chimney is 20 m high. Neglecting any plume rise, find the ground concentration 10 km downwind the refinery for the following conditions: (i) wind of 2.5 m/s at night on a clear night; (ii) wind of 8 m/s on a sunny day with strong sunlight. Calculate the horizontal width of the plume at this point. If the atmospheric temperature is 280K, is the assumption of negligible plume rise valid? Solution This is a straightforward application of Eq. (4.19). The mass flow rate of the pollutant at the source is Q=(100x106) x (/287/300) x (0.22/4) x 10 kg/s = 3.65 x105 kg/s. The dispersion coefficients for the given conditions at x=10 km are given from Table 4.1 and Fig. 4.4 as: (i) category F, so y=280 m and z=50 m; (ii) Category C, so y=800 m and z=530 m. Therefore, from Eq. (4.19), we obtain that the ground concentration is (i) 3.1x1010 kg/m3 and (ii) 3.4x1012 kg/m3. The horizontal width of the plume at that location can be defined as twice the y-coordinate where the value of the pollutant is 1/e of its value along the centre. This occurs at , which gives the overall width of the plume as 792 m for case (i) and 2263 m for case (ii). To calculate the plume rise, we use the empirical formula Eq. (4.22) which gives that hpr=1.34 m for case (i) and 0.42 m for case (ii). Hence, we were justified to neglect it, more so for the strong wind case (ii). 22 2 yy  This example demonstrates the major alterations to the plume dispersion process expected from a switch from stable to unstable stratification conditions. Example 4.5 For the plume of Example 4.4, find the point immediately below the axis of the plume on the ground (i.e. y=0, z=0) that has the highest pollutant concentration. Assume that and my Ax n z Bx . Solution Differentiating Eq. (4.19) with respect to x and using the given variation of the dispersion coefficients with x, we get that 0/)0,0,( dxxd when n mz z 22 1 2 2 2 0   . It is a reasonable approximation to assume that m/n is not substantially different from unity. Hence, the maximum will occur when 2/0zz  . Putting numbers, the maximum concentration will occur when z is 14.1 m. So reading from Fig. 4.4, for case (i), this will occur at a downwind distance of 1.5 km, while for case (ii), the maximum ground concentration will be at about 0.25 km. Note the very different behaviours between the two atmospheres, which shows the different possibilities between pollution dispersion during day and night at the same geographic location. 51 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment Appendix: An alternative derivation of the Gaussian plume equation For simplicity of presentation, we will consider only one-dimensional spreading. Assume that Q kg/s are emitted over a line of length H parallel to the y-axis. There are no solid boundaries. The governing equation for the mean concentration then becomes 2 2 z K x U       (A4.1) which is identical to the unsteady heat equation if we recognize that t=x/U. Our initial and boundary conditions will determine the solution. Many solutions of Eq. (A4.1) are given in advanced heat transfer books. An alternative way to show the rationale behind the solution we used (Eq. 4.20), is to seek self-similar solutions. By “self-similar” we mean solutions that depend on one independent variable only, which is a combination of the two dependent variables t and z. Most slender turbulent flows are self-similar for the velocity and the scalar fields. Let us say that we seek a function f(), such that )(0  f (A4.2)  /z (A4.3) with and  being a function of x only. The quantity 0 0 is a characteristic scale of the concentration, while  is a characteristic width of the plume, both as yet undetermined. We will use the following transformation rules: x f x f            1      z and from Eq. (A4.3), xx         Using the above and substituting Eqs. (A4.2) and (A4.3) in (A4.1) gives:                    2 2 2 000    d fdK d df dx dU f dx d U 00 0 2 2 2                 f dx d K U d df dx d K U d fd     (A4.4) If we want self-similar solutions, the terms in brackets in Eq. (A4.4) must be independent of x, i.e. constants. Let us consider the term dx d K U  first. Without loss of generality, we could say that this constant is unity (if it were anything else,  could be altered so that it is). This results in: 52 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment 5. Turbulent reacting flows and stochastic simulations 5.1 Fast and slow chemistry Turbulence time and length scales The turbulence has a wide range of scales. If the rms velocity fluctuation is u and the large eddy length scale Lturb, the large eddy turnover time, Tturb, is given by  u L T turbturb   (5.1) and the eddies at the Kolmogorov level have a time and length scale of (5.2) 2/1Re  tturbK TT (5.3) 4/3Re  tturbK L where  / Re turbt Lu  (5.4) is the Reynolds number of the turbulence. For a typical location low in the atmospheric boundary layer, we may have that u= 1 m/s, Lturb = 100 m, and for air =1.5x105 m2/s. These values give Tturb = 100 s, Ret = 6.7x10 6, and hence TK = 0.04 s and K = 0.8 mm. So there is a very wide range of time and length scales in typical flows in the environment. Chemical timescales – the Damköhler number From the presentation of chemistry in Chapter 2, it is evident that the rate of the reaction is determined by the amount of reactant and by the reaction rate constant k. The governing equation for a scalar (Eq. 3.9) has terms due to the fluid mechanics (advection, turbulent diffusion, molecular diffusion) and due to the chemistry. Therefore the net rate of change of the scalar will depend on both phenomena. Now, added to this that the various terms in Eq. (3.9) have different scaling even without chemistry (large-eddy scaling for the mean and turbulent transport, small-eddy scaling for the molecular diffusion), we see that the relative effect of the turbulence on the reaction will depend on the ratio of turbulent to chemical timescales, called the Damköhler number. We usually define two such ratios, one using the large-eddy timescale and another using Kolmogorov timescale: 2/1121 Re ,  t chem K chem turb Da T T Da T T Da (5.5) We do this because the chemical reaction occurs at the small scales, where the reactants mix at a molecular lever, but we are also interested to know how the reactants get together in the first place, and this is determined by the large scales. These complications are typical of turbulent flows. It is not just a coincidence that turbulent reacting flows are considered one of the most difficult problems in the whole of turbulence and fluid mechanics fields! 55 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment A way to define the chemistry timescale is illustrated by the following example. Take the reaction NO + O3  NO2 + O2. The consumption rate of [NO] from this reaction is d[NO]/dt =  k[NO][O3]. Treating [O3] to be constant, we obtain d[NO]/dt = A [NO], which gives an exponential decay with timescale defined as 1/A. Typical concentrations of NO and O3 in a typical early morning are 3x109 and 1.5x109 kmol/m3 respectively and the rate constant is k=1x10 m3/kmol/s. Hence, for this reaction and for these concentrations, Tchem=66.7 s. This is a fast reaction for atmospheric chemistry; others are many orders of magnitude slower. The definition of chemistry timescales becomes more complicated if there are many elementary reactions that consume or produce a particular species, as is almost always the case. This makes matters very difficult, as far as theoretical predictions go. Simplifications Fast chemistry (Da2 1): If the chemistry is considered to be very fast compared to all fluid mechanics scales, then the rate of chemical reaction follows the rate at which reactants are being mixed. This is the basis of many models, such as the “eddy break up”, the “mixed is burned” and the “eddy dissipation concept” (Peters, 2000). All these have been used extensively in the combustion literature, where indeed the chemistry is usually fast compared to mixing times, e.g. inside engines. The situation is less clear concerning environmental flows. There, the large eddy time scales are longer than in engines, but the Ret is higher, and hence the Kolmogorov time scale may not be very different. The chemistry, on the other hand, is much slower in the atmosphere than in a flame. This makes the fast chemistry limit invalid for atmospheric pollutants. Hence, the well-proven techniques that are embodied in commercial CFD codes, usually tuned to combustion or other chemical engineering problems, may not be always applicable to environmental reacting flows. Slow chemistry (Da1 1): At the other extreme, if the mixing times are much faster than the chemistry, we may consider that the reactants are well-mixed, i.e. that the pollutants are uniformly distributed in space, at least over a scale comparable to Lturb. This means that the degree of scalar fluctuation is small, which enables us to neglect the correlation 21  . This then allows the straightforward application of Eq. (3.19) without any need to provide closure for the turbulence effects on the reaction, which is calculated simply using the mean concentrations (an assumption sometimes called “laminar chemistry”). This assumption is very often introduced in air quality models that are currently used for monitoring and certification purposes. However, as these models progressively look into finer and finer length and time scales (e.g. the NO emitted from cars in street canyons), the assumption of well-mixed reactants becomes invalid. See Example 5.2 for more discussion on this. Intermediate (Da2  1 Da1): This is the most difficult, and the most relevant, situation, which is currently under intensive research. No simple theory exists for this regime. One very useful and insightful technique is to calculate the joint pdf of the scalars by numerical means, and this is introduced next. 56 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment 5.2 Monte Carlo simulations Fundamentals of Monte Carlo techniques In Monte Carlo techniques, a large number of events are sampled from the possible sample space and their evolution is tracked. Averages formed over many such “trajectories” give the behaviour, in a statistical sense, of the whole system. Other, perhaps more accurate, names for this technique include “stochastic modelling” or “stochastic simulation”. Monte Carlo techniques are used for an immense variety of topics (e.g. statistical mechanics, molecular dynamics, chemistry, turbulence, stock market risk analysis). The following examples aim to clarify the technique. Example 1 of Monte Carlo simulation Assume we have a die. Each time we roll it, we note the result (1 to 6) which we call x. Then, we use this number to initiate the following sequence of events: (i) calculate y(x) = (x/6)2 + sin(x/6) + x exp(x/6) ; (ii) if the result is less than 0.5, it is “YES” ; (iii) if the result is greater than 0.5, it is “NO”. We want to calculate the probability of “YES”. Clearly, an analytical approach is extremely difficult (but not impossible for this particular example). The Monte Carlo approach is the following: (i) sample N random numbers from a uniform distribution between 1 and 6; (ii) calculate y(x) for each of these N samples; (iii) for each of the resulting y(x), determine if it gives “YES” or “NO”; (iv) count the number NYES of events that were “YES”. The probability of getting “YES” is then simply NYES/N. All that is needed is a computer and a random number generator (and there are plenty of library codes in Fortran, C, or in Matlab for that). Example 2 of Monte Carlo simulation Assume you have £100.00 that you want to invest in 10 different stock market companies, S1, S2, …, S10, each of which today costs A=£1.00. History suggests that the probability of S1 rising tomorrow by a fraction f1 (i.e. the value tomorrow is equal to (1+ f1)A) is given by a normal distribution with mean 1=0.01 and rms 1=0.1. So, on average, the price of S1 tomorrow will be £1.01, but it could be anywhere between (1.013·0.1) and (1.01+3·0.1) (to within 99%). Similar historical data are assumed to apply to all our equities. There could also be various restrictions, for example: if it happens that S1 becomes greater than £1.50, then S3 will fall by 10%, etc. We want to calculate the spread in the possible worth of our investment portfolio after one month. A Monte Carlo simulation would be to simulate a very large number of possible “scenaria” and then evaluate the worth of the portfolio at the end of the month. The steps would be: (i) pick randomly a value for each fi (i.e. for each stock), each according to its own distribution (assumed known from historical data). Hence you get the new price of each stock Si. (ii) Enforce the various restrictions (e.g. change the Si’s accordingly). (iii) Continue for the next day. After doing the simulation for 30 days, you get the value of your investment for this particular realization. (iv) Repeat steps (i) to (iii) for N (some thousands) of these scenaria. At the end of the simulation, we have N possibilities, each with a different total worth of our investment. These define the pdf of our investment value, which helps answer questions like “what is the probability that we will lose 10% of our money” or “what is the probability that we will gain 20%” etc. Warning: this is a hypothetical example: DON’T TRY IT – no responsibility is borne by the lecturer if you lose your money! 57 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment z-z0 Z (a) y Y (b) Figure 5.1 The meandering plume. (a) View from the side. (b) View from above. Figure 5.2 The mean concentration across a meandering plume for different values of YRMS/y. Comments If we wanted to approach this problem analytically, the mean concentration at y would be given by:   dYYPYyUHQdYYPYyy yy )(2/)(exp)2/()(),()( 22        (5.8) 60 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment This is quite difficult to evaluate in closed form and if at all, it can be done only for simple forms of the P(Y), the pdf of Y. It can be shown that for a Gaussian P(Y), )( y is given by an expression identical to Eq. (4.20), but with an effective dispersion coefficient 222 , RMSymeandy Y  (5.9) This is typical of random, uncorrelated motions: the overall variance due to the action of two phenomena (meandering, smaller-scale turbulent diffusion) is the sum of the variances due to each phenomenon alone. (This concept is also used in experimental error estimation.) The basic assumption behind the flapping model presented previously is that the turbulent diffusion motions resulting in Eq. (4.20) are uncorrelated from the large-scale motions causing the meandering. However, in real atmospheric turbulence, the scales causing mixing of the plume are not fully independent from the motions causing the meandering. Hence, the discussion above should be seen as only an approximation to a real plume exposed to sudden, strong, cross winds. Equation (5.8) cannot be easily evaluated analytically for other shapes of the pdf of Y. However, the stochastic simulation does not have this restriction and it also gives access to the full pdf of , which is very useful for predicting the fluctuations of pollutant at a particular point. 5.4 Worked examples Example 5.1 On a particular day, the air above a city is completely still. A box model (Section 4.1) is here extended to account for small-scale turbulent fluctuations in the concentrations that are a function of time. Assume that there are no emissions and that the pollutant is inert with homogeneous mean concentrations in all directions, but that there are significant fluctuations at t=0. Hence, simplify and solve Eq. (3.16). The turbulence lengthscale is Lturb=50 m and the turbulence velocity is u = 0.1 m/s. Calculate the rms after 5 min as a fraction of the initial rms 0.  Solution When there are no spatial inhomogeneities, Eq. (3.12) becomes consttzyx dt d w x D xx u t j T jj j                    ),,,(0   , as expected. The equation for the variance (Eq. 3.16, using Eq. 3.17) becomes: 2 2 ' 2  turbL u dt d  . Integrating this equation gives: , with  turbTt /2exp202   '/ uLT turbturb  . This is an exponential decay to zero from the initial value of 0. This shows that the inert scalar fluctuations will have decayed to 1/e their initial value over one half turbulence timescale and that they will have essentially disappeared after three turbulence times. For the conditions given, Tturb=500 s, and hence  /0=0.3. Note the quick disappearance of the fluctuations in atmospheric turbulence, e.g. of the order of minutes. This is consistent with everyday experience. 61 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment Example 5.2 Repeat Ex. 5.1, but now assume (a) that the pollutant is photolysed at a rate –Jc, where J is a constant and c the pollutant concentration; (b) that the pollutant is destroyed at a rate –Jc2. Solution In this Example, the mean concentration is changing with time due to the chemistry and we need to include the fluctuations in the chemical source term. (a) The governing equation becomes:  JtcccJJc dt cd w dt cd  exp)0( . The corresponding equation for the reacting scalar variance is now        J T ccJ T wc L u dt d turbturbturb 1 22 2 2 ' 2 222 2   . The solution is easily obtained as: . So now the decay of the fluctuations depends on the turbulence and   )2exp(/2exp202 JtTt turb   on the chemistry. If the initial mean concentration is C, an appropriate chemical timescale is C/J. We can then define a Damköhler number as Da=JTturb/C. If C/J>>Tturb (i.e. Da<<1), we have intense turbulence and the fluctuations will have decayed to zero well before the reaction is felt. On the other hand, if C/J<<Tturb, the fluctuations decay because the reactant is consumed ( ) according to the chemical timescale and the turbulence plays little role. This Examples clearly demonstrates the Damköhler number criteria in Section 5.1. 0c (b) In this case, we have non-linear chemistry and the system of equations becomes: 22 JcJ dt cd ccJ dt cd w dt cd   cccJ Tdt d wc L u dt d turbturb  22 2 2 ' 2 2 2 2 2   32 2 4 2 cJcJ Tdt d turb         . (   32222 2)2( cJccJwccccccJwwccJwwJcw   ). It is evident that we cannot proceed unless we supply a model for the third-order correlation. This example shows how difficult the combined turbulence-reaction system can become when the chemistry is non-linear. However, it also demonstrates the very useful conclusion that for fast turbulence relative to the reaction ( ), we can safely assume that the reactants are well- mixed ( 0Da 0 ). Virtually all Air Quality Models in practice are based on this assumption, which is used without much thought even if Da is not very small due to the large simplification it offers. 62 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment Stability categories and empirical data on dispersion coefficients: A to F refer progressively from the very unstable to the very stable conditions, with D the neutral. Day Solar intensity Night Wind speed at 10 m high (m/s) Strong Moderate Slight Overcast Clear 0-2 A A-B B - - 2-3 A-B B C E F 3-5 B B-C C D E 5-6 C C-D D D D > 6 C D D D D Horizontal 65 4A8: Environmental Fluid Mechanics Pollution dispersion in the environment 66 Vertical