COO-1407-39 THE UNIVERSITY OF MICHIGAN COLLEGE OF ENGINEERING Department of Meteorology and Oceanography Progress Report No. 8 RAIN SCAVENGING STUDIES A. Nelson Dingle ORA Project 068670 under contract with: U. S. ATOMIC ENERGY COMMISSION CONTRACT NO. AT(11-1)-1407 ARGONNE, ILLINOIS Adminiatered through: OFFICE OF RESEARCH ADMINISTRATION ANN ARBOR January 1972

~i1^1b ", P- X

TABLE OF CONTENTS Page LIST OF TABLES AND FIGURES iv I. INTRODUCTION 1 II. EMPIRICAL STUDY OF THE FIELD DATA OF 1 JUNE 1970 7 A. Detailed Description of Deposition Pattern 8 1. The Axis of Maxima 8 2. The Mass Deposition per Meter 8 3. Forking of the Axis of Maxima 10 4. Inference of Cyclonic Rotation 10 5. Relation of Deposition to Concentration Patterns 13 6. Separated Maxima of Deposition 30 B. Extraction of Quantitative Information 34 1. Axis of Deposition maxima 34 2. Estimation of Tracer "Plume" Concentrations 38 3. Longitudinal Profile of Deposition Mass 44 4. Inference of Deposition Coefficients 46 III. THEORETICAL MODELING 51 A. A Brief Review of Entrainment and Numerical Experiments of Cumulus Clouds 52 1. Entrainment models for Cumulus Clouds 52 2. Numerical Experiments 57 B. On the Vertical Distribution of Contaminant The Model 64 1. Introduction 64 2. The Model 64 3. The Updraft Radius 65 4. The Updraft Trajectory 67 5. The Contaminant Distribution 70 6. Computations and Results 73 7. Conclusion 86 C. Cloud and Rain Formation 86 1. Introduction 86 2. Assumptions 86 3. Equation of Motion 87 4. The Entrainment 87 5. Equation for Water Substance 90 6. Thermodynamics 96 7. The Environment 98 8. Rainfall Rate 99 PUBLICATIONS 102 PERSONNEL 103 iii

LIST OF TABLES AND FIGURES (Tables and Figures from Part I. and II.) PAGE TABLE 1. Timing data for scavenging coefficient estimates 25 TABLE 2. The relative variation of r along streamer of Tracer-In. 41 FIGURE 1-a.. Distribution of tracer indium concentration, ng 9Z-. 3 1-b. Pattern of tracer indium deposition, ng m-2. (see figure 1-c for values). 4 I-c. Distribution of tracer indium deposition, ng m-2. 5 2. Mass of indium deposited per kilometer of di tance along the deposition streak axis. 9 3. Rawinsonde for Peoria, Illinois, 1900 CDT, 1 June 1970. 12 4. Rawinsonde for Salem, Illinois, 1900 CDT, 1 June 1970. 14 5-a. Radar echo patterns, 1639-1650. Also shown in the indium flare burn (1639-44). 15 5-b. Radar echo patterns, 1651-1745. 16 6. Indium concentration and rainfall rate as functions of time at Michigan Base Station (M). 19 7. Indium concentration and rainfall rate as functions of time at Automatic Station no. 1 (G1). 20 8. Indium concentration and rainfall rate as functions of time at Automatic Station no. 2 (G'2). 21 9. Indium concentration and rainfall rate as functions of time at Automatic Station no. 3 (G'3). 22 10. Indium concentration and rainfall rate as functions of time at Mobile Station no. 1 (Mi1). 23 11. Indium concentration and rainfall rate as functions of time at Mobile Station no. 2 (M2-). 24 12. Map showing estimated travel of tracer plume in relation to sequential stations. Also indicated are the times of approach of leading edge (ti) and of passing of trailing edge (tf) estimated by various schemes. 29 iv

LIST OF TABLES AND FIGURES (Cont.) FIGURE PAGE 13-a. Tracks of the center of mass of radar echo cells A, B, and C. (see also figure 5-a). 32 13-b. Tracks of the center of mass of radar echo cell D. (see also figure 5-b) 33 15. Lateral deposition profiles at several distances along the deposition streak. 37 16. Trace deposition along streak axis. 38 17. Comparison of values of X and. D along the deposition streak axis. 43 (Figures from Part III.) 1. Verticle wind profile. 74 2. Mean draft radius. 75 3. Horizontal wind speed. 77 4. Draft Shape. 78 5. Draft Shape. 79 6. Draft Shape. 80 7. Vertical profile of contaminant in a disc of thickness AZ = 100 m. 82 8. Vertical profile of mean mixing ratio in a disc of thickness AZ = 100 m. 83 9. Radial profiles of mean mixing ratio. 84 10. percent of mass of contaminant in a cylindrical shell of radial width 500 m and of vertical thickness 100 m. 85 V

I. INTRODUCTION In approaching the present report, it appears useful to review briefly the background of our experiments with indium tracer in convective storm systems. An early objective was stated in somewhat over-simple terms as follows: to identify the beginning and ending points of tracer particle trajectories in connection with convective storm activity. The point of this objective was that, if achieved, it should provide a basis for evaluation of the intervening processes (rain formation, circulation and scavenging) within a somewhat more limited domain than previous observational data could supply. The jungle represented by "rain formation, circulation, and scavenging processes" is exactly that which we are exploring. By first examining our experimental data critically in detail we have assured the most reliable possible definition of the limits of our speculative domain. The map of indium deposition of 1 June 1970 (Figure lb) represents this definition for our most informative experiment. The fact that this does not contain all the ending points of the indium particle trajectories is in itself informative, though not surprising. Obviously the limits of our sampling network do not comprehensively contain the detectable deposition streak. To wring the maximum of information from these data, we have undertaken to study in some detail both the deposition pattern itself and the processes that must intervene between the initial emission and the deposition of the tracer material. It is this study which will finally define with greatest resolution the physical processing of the tracer aerosol from which "models" of the experimental storm may be constructed. 1

From this point, more general models may be derived and more definitive experiments designed. Although these remarks suggest a temporal priority for the empirical phase of our research, the best use of our manpower and talent has indicated that simultaneous pursuit of our theoretical modeling objective is in order. Accordingly, the main body of the present report is divided into Empirical and Theoretical sections. The third section deals with ancillary studies and sub-projocts. 2

CODE: N- ISWS STATION N N - UM STATION 2.4 16 (N)-SEE TEXT -30 KM -7 i 3 9 T9 T T - -.1 -- - - --! 3:9; OT 16:3 34.8 70 1J.-, T'" "'"" I I T T T - T i I f f'i.. ^.. 7:0 3 2;.51. 1 37.:8 -3: 32.2 A T 3:24- 44 -- -- 4- - - i 15 K.. r T I. 229. 12.6 - T T 45K5 763.7 T T 17R0 T m T 5 /, 42.9 I 9.f * I T 310 615 (T7Z) 37:9 5.4 I\ 148 229. i i f, 2 4 6 8.0 KILOMETERS -0 KM f T f T BURN TRACK 1900 FT MSL TRACER IN CONCENTRATION -NG/1 (TOTAL MINUS BACKGROUND OF 4.0 NG/L) I JUNE 1970 FIGURE i-a: Distribution of tracer indium concentration ng p-'. 3

\TRACER IN DEPOSITIONN 3015 KM ) /' 8\ 10 \\TRACER IN DEPOSITION-N/1 (TOTAL MINUS BACKGROUND OF 4D0 NyL) -0KM 9 X/O^1 I JUNE 1970 I0 I7 1 X0 1g 1 ISOPLETHS IN POWERS OF 2 2s-32 2"'2048 BURN TRACK 1900 FT MSL 0 2 4 6 8 10 I 1 JUE17 I00~~ ~KILOMETERS FIGURE l-b: Pattern of tracer indium concentration, ng m-2. k~~~8-,3 24

75 t CODE: N- sws STATION N N - UM STATION 30 62 (N)-SEE TEXT -30 KM. 240;51; _o,00 T I' T T _ _ I I 490 10 220. 5 IJ T T 1500. I;9;0.3i 2 - | 7____-220-570 ----,20 T. 7o o T o 15 KM 5 - - T 4 10 ~620 30 * TOT o -ToT 00 100 0.T T ~i. T ( T T /*900 * 250,.. I. T. 48 1050 (90) 50.i T I - T T 460 io.. I 830 450 350 300 14bO i00 i i - 3/00 ____ T T T T o, 450 KILOMETERS I JUNE 1970 FIGURE l-c: Distribution of tracer indium deposition, ng mn2 5

II. EMPIRICAL STUDY OF THE FIELD DATA OF 1 JUNE 1970 7

A. DETAILED DESCRIPTION OF DEPOSITION PATTERN THE AXIS OF MAXIMA The observed deposition pattern (Figure l-b) may be described more or less quantitatively in terms of the decrease of deposition per unit area, D, along the axis of maxima. This decrease is approximately exponential, and is reasonably estimated by D =D e (1) where D is the deposition in ng m-, D is that at the edge of the network, x is the axial distance in km, and c = 0.069315 [=-.l(-ln2)]. (i.e., in 10 km, D decreases by a factor of 2.) THE MASS DEPOSITION PER METER An additional indicator for the rate of decrease of the deposition with distance from the initial maximum is found in the total mass of tracer deposited within the deposition streak. To evaluate this, the mass M of tracer deposited in a strip 1 m wide taken normal to the axis of maxima was computed: dM Dd M= d= Ddy (2) These values are plotted in Figure 2. Clearly some additional information is found here, and new interpretations are needed. We derive the following from these data: d for 0 = x < 6 km, a-(lnM) = -0.105 d for 6 < x <14.5km, d-(lnM) = 0 8

20108- 0 8 o 4- o - 0 10 20 30 X (KM) MASS OF INDIUM DEPOSITED ALONG THE DEPOSITION STREAK AXIS I JUNE 1970 MX = TOTAL MASS FROM 0 TO X FIGURE 2: Mass of indium deposited per kilometer of distance along the deposition streak axis. 9

for l4,5t x < 30 km, a(lnM) =-0.051 Although the plateau from 6 tol4.5km appears troublesome, the possibility that the different deposition rates are affected by different scavenging mechanisms emerges here, and further exploration is called for. FORKING OF THE AXIS OF MAXIMA Figure l shows a splitting of the axis of maxima beginning near 15 km and resulting in a two-pronged principal pattern farther downstream, with additional details that appear to be attributable to the cellular character of the precipitation field. The two-tongued pattern has been observed in other field experiments (Chem-Met Workshop, 1971). The parallelism of the two nearly equal axes of maxima suggests parallel-moving cellular systems. The dynamical scales and details remain to be investigated, probably most effectively with the aid of radar echo analysis. We have commented previously (Dingle, 1970) upon the sharpness of the right side of the tracer concentration pattern. In the present instance, we have made detailed studies of the anomalous-looking stations, leading to corrections of the earlier presentation. The deposition pattern in its final form (Figure l-b) again shows the sharp gradient east of the axis of maxima, and a more diffuse gradient on the west side. INFERENCE OF CYCLONIC ROTATION Subjective isopleth analysis can only be that; however, using the greatest care and applying the technique so as to express the data as accurately and completely as possible, 10

one cannot avoid the implication of cyclonic curvature of the tracer trajectories. This also was noted, though less positively, in a previous experiment (Dingle, et. al., 1969). No comparable implication of anticyclonic trajectories is observed. An inferred cyclonic rotary tendency of the updrafts in individual convective cells is further supported by the available wind data (in the PIA and SLO soundings, Figure 3 and 4). At the time of the soundings, SLO was east of the front, and should represent best the flow and stratification in the lower atmosphere over the sampling network; PIA was west of the front, and probably represents the upper level flow patterns best. The low level flow was very nearly from South; a slight tendency for SSE winds at the surface was observed at a few stations. Veering with height is indicated up to ~ 3 km at SLO, and very slight backing from 3 to 9 km, giving a 55 kt S wind at 9 km (top of the wind sounding). At PIA the veering effect was only slight but continuous to 9 km where the 45 kt wind blew from 9 190~. These winds cannot be combined in any way to produce transport of the tracer particles from the burn track off the western edge of our sampling network; yet only one sample along that edge was below the background criterion for definite tracer deposition. We therefore are forced to the conclusion that the wind soundings do not at all indicate the winds-inside the convective rain-generating cells. Further, we are forced to the inference that these cells did contain major cyclonically-turning currents which were capable of transporting substantial amounts of tracer material at least 8 km west (left) of the dominant general wind trajectory. The radar data have been traced from 35-mm b-w film copied from the original and supplied to us by the Illinois 11

P U IO(Mb) (Niv) / / / / / (.,.b)/300, 30 -40 I 25' 400 o Td 2 /0 -30~6 —20 6Tc 500 51 15 600 - -20 N. /\ \_ /700 3-10,,00 2-'5 -10 0 10 20 0-L0 T(oc) PEORIA I JUNE 1970 1900 CDT FIGURE 3: Rawinsonde for Peoria, Illinois, 1900 CDT, 1 June 1970. 12

State Water Survey. Although a stepped gain procedure was followed, we find it virtually impossible to put this technique to use in the 1 June 1970 analysis. Figure 5 shows some of the radar echo patterns in their time sequence and in relation to the network stations. An immediate caution in interpreting the echo patterns is in order: the echoes represent only the field of water particles large enough to return the radar signal, and they are located at altitudes of 300 to 500 m. They are therefore closely indicative of rain received at the ground, but their relation to updraft and other circulational details is not clear. The largest echoes have dimensions of the order of 3 km by 5 to 7 km, the long axis tending to be nearly N-S oriented. Recalling that these are areas of falling rain, that is, generated at a higher level, it is likely that the most vigorous updraft currents are peripheral to the echoes. The echoes themselves, therefore, can neither confirm nor refute the hypothetical cyclonic circulation, but the initial and endpoints of the tracer-In trajectories appear definitely to require it. RELATION OF DEPOSITION TO CONCENTRATION PATTERNS It is necessary to clarify the relation of observed values of tracer concentration to the derived deposition values. We are involved here with considerations of the field conditions (temporal relation of tracer passage to the occurrence of rainfall, for example) and with problems of analysis in the laboratory. Specifically, the neutronactivation analysis is a measure of the total indium contained in a sample of measured volume taken from the rain catch at each station. It is therefore a measure of the concen13

p - 10o (Mb) (W/s) 300.30 9 -40' /T /!(/ w 2 ('Or. /400 //(-20 0 -30 6 35-15 600 -20 700 - 1 / ______ /_ _______1_ _______0_ _______ / ~^^ _____ ^/ 1000 - 0 1 20 T(oc) SALEM I JUNE 1970 1900 CDT FIGUIRE 4: Rawinsonde for Salem, Illinois, 1900 CDT, I June 1970. "' - 5 //~~~~~~~~~~~~ \ ~ ~ ~ ~ ~ ~ ~~000 0 ~ ~~ ~ ~~~~~~~~O 102 T(7/4-c) -=~ 0 ~ ~ SLE 1 JUNE 1970~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~: 7' r1900"CDT FIGUE4 aisnefrSlm lios 90CT 1 June 1970~~~~~~., /~ ~oo ~-,o~~1

1.. 0 10 KM 1639 CDT 1640 1641 1642 ~ ~ ~ ~ ~. 0~ ~ ~.. * 0. 1643 1644 1645 1646 B QO 0 0 1647 1648 1649 1650 INSTANTANEOUS RADAR ECHOES I JUNE 1970 FIGURE 5-a: Radar echo patterns, 1639-1650. Also shown in the indium flare burn (1639-44). 15

GCC N~ ~.,O.. B... * ^10 KM 1651 CDT 1655 1700 1705 ~ ~.s ~.. ~.'1o.I....~c...(7 ~ ~ 0D0 D,, 1730 1735 1740 1745 INSTANTANEOUS RADAR ECHOES I JUNE 1970 FIGURE 5-b: Radar echo patterns, 1651-1745.

tration, X, of the tracer in the whole water sample. This then is the basis upon which the background level is measured. Whereas in Oklahoma, we found a background level near 7 ng V"1 (Dingle, et al., 1969), numerous measurements in the Illinois network area indicate a background of 4 ng -1 Thus tracer indium is identified "with certainty" only in those samples for which the concentration exceeds 4 ng a', regardless of what the total deposition may have been. It is therefore necessary to delimit the identifiable tracer deposition streak by means of the 4 ng'"1 concentration isopleths. This gives a reasonably distinct border along the east side of the deposition pattern (Figure 1), but is outside of the network to the west. To evaluate tracer indium deposition, therefore, we have deducted 4 ng -l1 from the concentration determinations (Figure 1-a), observed the border given by the 4 ng -l1 isopleth, and computed deposition amounts for stations west of this border, using D = Q X (3) where D is deposition of tracer in ng m-2, Q is total rainfall in mm, and X is rainwater concentration of tracer indium in ng V'1. Because the whole-storm samplers simply accumulated rainwater until precipitation ceased, the tracer-tagged portion of the rain was generally diluted considerably. For purposes of future experimentation it should be noted that the size of sample taken for indium analysis should be proportional to the total rainfall at each station. Had we followed this procedure in place of the constant volume (1 9) sample procedure, we should probably have improved our estimation of the tracer deposition pattern. In any case, 17

the basic figure derived from the neutron activation analysis is a diluted average rainwater concentration which must be multiplied by the total rainfall collected at the station to give the deposition figure. We do not, therefore, have direct information on the concentration of tracer in the falling rain for the whole-storm samplers. We have such data only from the sequential sampling stations. Figures 6-11 give the results of the indium measurements for the sequential samples. These were acquired at 6 stations (see map, Figure 12) 1) Michigan Base Station at the northern edge of Clinton, designated M 2) Automatic Station #2, about 4 km west of M, desianated G*2 3) Automatic Station #1 about 2 km southwest of, designated G 1 4) Automatic Station #3, about 2 km east of M, designated G*3 5) Mobile Station #1 located at Rainguage station #104 about 12 km NNE of M, designated M1 6) Mobile Station #2 located near Rainguage Station #77 about 24 km NNE of M, designated M2. The normal background concentration level of 4 ng Z-lapplies also to these data. We do note, however, that the effect of dilution that is normally observed in sequential samples as the rainfall rate increases (Bleichrodt, et al., 1959; Dingle & Gatz, 1971; Georgii, 1969; and others) appears to justify an identification of tracer-In with lower concentration values in some of our sequences under relatively heavy rain. For example, at M (Figure 6) the first shower appears definitely to be labeled with tracer-In in its light-rain onset and decay periods whereas the two samples collected at the climax of the shower fall below the 4 ng 7-1 threshold. The deposition rate computed for this sequence shows a maxi8

INDIUM (ng/l) TRACE AMOUNTS i\ ^POLL~NPEN F a (grains/ml) 100: \^ /^^______10 I _ - B RADIOACTIVITY I 12 (pc/I) RAINFALL RATE (mm/hrn) 30 -O 3 - CDT JUNE 1, 1970 MICHIGAN BASE FIGURE 6: Indium concentration and rainfall rate as functions of time at Michigan Base Station (M).

100 INDIUM (ng/l) 10 180-,20RAINFALL RATE (mm/hr) 60-'0 -_\\ i2so 75 1730 I S 18;00 81 CDT JUNE.,1970 AARC -I FIGURE 7: Indium concentration and rainfall rate as functions of time at Automatic Station no. 1 (G*1). 20

100 INDIUM (ng/l -10 AT 2 03+ TRACE AMOUNTS 120- -. 80^ RAINFALL RATE (mm/ hr) L0O ":~~~. 13~ ~ * 17I0 1800 130 9I I I i 1350 1rOO 1730 1800 1830 1900 2200 CDT JUNE 1,1970 AARC -2 FIGURE 8: Indium concentration and rainfall rate as functions of time at Automatic Station no. 2 (G'2). 21

10 INDIUM (n//l) t6'~~~ 0-^ ^AMOUNTS 0-12 8 RAINFALL RATE (mm/hr) 0- I 1250 1730 \7 1i08 1815 1830 CDT JUNE 1,.1970 AARC-3 FIGURE 9: Indium concentration and rainfall rate as functions of time at Automatic Station no. 3 (G'3). 22

\^ TN CONC. \S~ I (n/9) tIoot RAINFALL RATE 40 20 1720 1740 1800 CDT,, JUNE MOBILE 1 FIGURE 10: Indium concentration and rainfall rate as functions of time at Mobile Station no. 1 (Mi*). 23

tO I. CONC. I - Ilos ~ ~ C~ II RAINARLL RATE I 40~' I t r 7~20 - lIzo?O0I7?o l8oo CDT, I JUNE MOBILE 2 FIGURE 11: Indium concentration and rainfall rate as functions of time at Mobile Station no. 2 (M2 ). 24

TABLE I: Chronology of Tracer-Tagged Samples at Sequential Sampling Stations - 1 June 1970. STATION G.1 Time Interval Rainfall Rate Concentration Sample Echo mean Speed mm/hr- Sample mean 1250- 17090.15 10.15 1715.6 1725 1715.61716.5 1 51.3 50.00 1716.5 171 7.0 81.8 26.36 1717.0 17171717.3 3 131.0 24.63 1717.31717.31717.6 1 149.0 26.04 177. 61717.9 1138.0 33.33 1717.9 177.2 - 146.0 17.79 1718.2 1718.2- 181.0 12.34 718.8 - 117.0 10.13 1718.81718.8- 1 27.9 9.67 1720.4 1720.4720.4 1 41.4 26.44 1721.4 1721.4- 3.8 31.80 1733.4 1733.4- 7 5 3.7 5.45 1744.2 25

TABLE I Continued STATION G*2 Time Interval Rainfall Rate Concentration Sample Echo mean Speed mm/hr1 Sample Mean 1704.85.4 4.10 1713.8 1710- 4 4 1713.836.6 36.91 1715.0 1715.0"1715.0- I 61.9 38.23 1715.8 171 6.4 111.5 90.91 1716.8171.8-; 61.2 73.19 1717.6 11735.0 1726 2.6 78.78 17351. 1749.0 3.4 67.69 1749. 0 1749- 18.4 18.70 1751.5 751.8 5 20.4 39.56 1753.8 175386 220.4 19.41 1756.0 1756.0-;30.7 3.56 1906.1STATION M 1715.0 1715:4 1714 4.85 1725.4- 2.85 1726.1 1 1726.1- 338 1727.1 1726.2-1 1730.5 1729 9.18 ~171530~.5-26

TABLE I Continued STATION G'3 Time Interval Rainfall Rate Concentration Sample Echo Mean Speed.mm/hr-' Sample Mean 1250- 1717 0.2 15.91 1734.3 1732 02 1734.313.6 4.1.2 (1.79) 1817.6 1817.61818.16 110.0 2.75 1818.1 STATION M1' 17251735.5 1.05 4.91 1735.5 1735.53.28 3.58 1746.7 1736 1746.71748.1 26.64 2.11 1748.1 1748.1749.2 1 29.30 2.05 1749.2 1749.2- 1749.7 1751 77.76 4.14 STATION M2~ 1718178 1755 0.70 8.39 1800 1800- I 1805 6.00 9.16 1805 1809 "1809 { 7.50 8.33 1809- LEAK 1810.8 I 1810.8181.8 1813 8.86 4.73 1814 181415.3 2.04 1816 181615.3 2.00 1818 1818- A I r, ~1818- 41.6 1.51 1818.8 27

imum for the most diluted (second) sample of 101 ng mI- hr-1. The fifth sample in the sequence, on the other hand, distinctly shows no tracer-In even though the rainfall rate during this sample period is very low, but the next two samples (second shower climax and decay, resp.) again indicate the presence of tracer-In. An independent estimate of the tracer-In position as a function of time may be constructed from the wind data. The earliest possible time of arrival of tracer-In at each station, ti, is inferred in first approximation by assuming straight line transport at an average speed of 35 kt (1.079 km min-) from the point of ignition of the flares. Isochrones of t. for 1650, 1700, and 1710 are shown in Figure 12. In particular, at M2., t. = 1717 CDT. The initial sample collected at this station (Figure 11) was collected from 1718 to 1800, so we are confident that no appreciable tracer-In was missed in the sampling. At G.l, however, the initial sample was collected from 1250 to 1716, and had an In concentration of 10.15 ng'1. But t. for this station is 1654, so we are confident that the entire tracer-In charge of this sample was acquired during not more than the last 22 min of this collection. A more precise timing is obtained from the radar record, as far as it goes. Figures 13 and 5 show the progress of identifiable radar echoes across the sampling network from 1639 to 1745 CDT. Regrettably, although our radio notes indicate that the radar was not shut down until 1815, the films supplied to us end at 1745. The most persistent (or persistently regenerated) storm cell (cell D, Figure 13 and 5) is associated with the occurrence of rainfall rates in excess of 100 mm hr-l, and gives a basis for estimating the rate 28

*M2 ~.1814 t =1710 CDT = 1747 *MI 1-1750 i =1700 t =1725 ti-1705' Go.2 -1734 t Gl;j m- 1715 =71710 tI-744 17- 17 40 t =1650 1t= 703 0 2 4 6 8 10 KILOMETERS \ = 1700 t' = 1718 TRACER PASSAGE TIME ESTIMATES: =13 t - FROM SEQUENTIAL SAMPLES t - ASSUMING 35 KNOT TRAVEL I \ tI = 1650 t- FROM RADAR ECHO MOTION |I \ ~= 1656 (0.5 KM MIN) = 1644 CDT ISOCHRONES OF PLUME ADVANCE I JUNE 1970 FIGURE 12: Map showing estimated travel of tracer plume in relation to sequential stations. Also indicated are the times of approach of leading edge (ti) and of passing of trailing edge (tf) estimated by various shcemes. 29

of progress of the tracer-In material along the eastern edge of the deposition pattern. Extrapolation of the recorded speed and direction carries this cell over M2 at 1804 to 1821 (t' and t;, Figure 12). The peculiar propagation pattern of this echo in the period 1725 to 1735 suggests strongly the shift from a mature and dissipating shower (1725) to a young and vigorous developing shower at the right rear flank of the older storm (1730-35) and then gains speed along the extension of the path of the earlier echoes. There is also a suggestion of accelerated movement associated with dissipation of the mature cells and deceleration associated with the intensification of the younger cells. These features suggest that this was not, therefore, a single steady-state severe storm (Browning, 1964), but was rather a self-propagating sequence of moderate-to-heavy thunderstorm showers (Newton, 1967). The temporal data are assembled for those sequential samples that appear to have contained tracer-In in Table 1. In most cases, the timing given by the radar echoes lies nicely within the time periods recorded for the tracertagged samples. In some cases the tracer-In appears to have lagged somewhat, continuing to arrive after echo passage. SEPARATED MAXIMA OF DEPOSITION Figures 1-a and-b show maxima that appear to depend upon isolated anomalously high values of In-content. Those along the eastern edge of the pattern appear to be wellexplained by the above discussion of the "self-propagating" strong shower cells that are shown by the radar analysis. The separated maxima along the western side of the pattern 30

evidently represent extreme penetrations of tracer-In material caught in the low-level influxes feeding the storm cells. We are dependent upon the ISWS rainguage records for the timing of showers in these cases, and from these and the requirement for transport from the burn track, we can infer the circulational features with some fidelity. 31

Mel N i G(3 G*2 A 1654, G-l 1654 TRACK C' G 1644 TRACK C 1654 CDT "TRACK A _ 1654 CDT Jst30 2 4 6 8 10 KILOMETERS 163 I RADAR ECHO TRACKS TRACK B-^ I JUNE 1970 I POINTS AND TIME MARKS INDICATE i APPROXIMATE CENTER OF CELL MASS I (^ /'I DOTS ARE AT ONE MINUTE INTERVALS 1639 CDT FIGURE 13-a: Tracks of the center of mass of radar echo cells A, B, and C. (see also figure 5-a). 32

I 351 MI )| / TRACK D2 30, 5 20 25 I~/ 30 G.3 IG 2 I Gel 15 20 10o I TRACK DI 15 TRACK D X / 1710 CDT 05 05 0 2 4 6 8 10 I -- -- 1'-' X --- - KILOMETERS 1700 CDT 57 i, RADAR ECHO TRACKS I I I JUNE 1970 1655 I POINTS AND TIME MARKS INDICATE APPROXIMATE CENTER OF CELL MASS 0~I FIGURE 13-b: Tracks of the center of mass of radar echo cell D. (see also figure 5-b) 33

B. EXTRACTION OF QUANTITATIVE INFORMATION AXIS OF DEPOSITION MAXIMA The axis of maxima of In-deposition, and its branches (Figure l-b) were first sketched in subjectively following the implications of the isopleth analysis. Cross-sections were constructed at various distances along the main streak axis for the purpose of inspecting points of special interest such as a) secondary maxima of deposition along the axis, b) the development of the forked pattern, c) the best location of the axis (Figure 15). The profile of the deposition maxima along the principal axis is shown in Figure 16. A first approximate expression for this curve is taken from the straight line: log 1D =log D- cx, c = 0.09 km2 x 2 0 ^ log D = 11.2 2 o = 11.2 - 0.09 x (4a) or A~ /% 1 Dx D e c, c1 =.0624 km- (4b) x o A more complete characterization of the data includes the pseudo-periodic variation that appears to change in amplitude gradually downstream: A ^ ~~~~ A^~~~ 27rx log D =log D + 0.5 Ax cos, = 12.4 km 2 X 2 X X (5) A =A + bx, A =0.6, b= 0.025 km-1 x o o 34

Conversion to the Naperian base gives 2~x In D = in D + In 2(0.5{A +bx} cos 2 - cx) x o 0 12.4 or D = D exp[(0.2079 +.00866x) cosL24 - ex] (6) X 0 1 2 ~ 4 The points indicated by x marks are computed from this equation. The value of equation (6) lies mainly in such relationships as may be developed between this expression and the physical processes that produce the deposition pattern. The linear decline of D, for example, may be most indicative of a mean scavenging effect operating upon a steadily diluting (diffusing) cloud of tracer-In. The uneven diffusion and the uneven scavenging to be expected in a field of rain showers offer a physical basis for the oscillatory behavior observed. There is, perhaps small reason to expect the constant value of 12.4 km for the wavelength of the oscillations. Considering that the radar echo mean speed is about 0.5 km minl, the maxima are about 25 min apart in time. This is a reasonable period to associate with the cycling of tracer plume material through a convective shower system. Because of the logarithmic scale, the increase of A x with distance does not actually indicate an increase of the amplitude of the deposition, Dx, with distance. The halfamplitude, a/2, of D may be studied as a function of x by noting that - = Do exp(0.2079 -.0537x) - D exp(-0.0624x) 35

dx (2) 57 0.2079-.0537x -0.0624x = -.0537 e +.0624 e D 0 Thus it is clear that d/dx (a/2) < 0 because the magnitude of the first term on the right is always greater than that of the second term, and the first term decreases in magnitude more slowly than the second term as x is increased. Hence the amplitude, a, of D decreases slowly with distance downwind, showing the effects of diffusion and depletion of the airborne tracer. The initial high value of deposition may most readily be attributed to a combined effect of dry fallout and impact collection by rain of the large aprticle component of the tracer-In plume. It is physically reasonable that this stage of the deposition process should decay rapidly to be followed by the deposition of the large-nucleus component. The second crest of D is found at a suitable distance x and time (12.4 km, 25 min) downstream to reflect an increasing contribution of nucleation scavenging, first of the largest nuclei, and later of smaller and more deeply penetrating nuclei. The third crest follows the sharp drop that is associated with the splitting of the axis of maxima. This suggests a broader dispersion and more rapid dilution of the airborne tracer cloud. The third principal scavenging mechanism, that of the diffusive attachment of the smallest particles to cloud- and rain-drops, operates as a relatively steady cumulative process throughout the period of association of the tracer and the cloud. It cannot be discerned separately by the inspection of the data. 36

12 10 -J 6-37 10 1 15 7 0 KM 3 5 7? 0 13 14,5 17 0 5 KM _~A^~. ~LATERAL PROFILES ALONG DEPOSITION STREAK FIGURE 15: Lateral deposition profiles at several distances along the deposition streak.

.~~~~2r 8p 72Dy-> B@-Ce 0+o.5, /cos C(2O f.) _A A o= IX; Aoo =., b=.o25 km'' A A A I1.2 _ -, D,,x =,n,2 o x 3, P,, =2 - ^o2 O - 0.09 x; D = (2) o 5 20 5s X(krn) TRACER DEPOSITION ALONG STREAK AXIS I JUNE 1970 FIGURE 16: Trace deposition along streak axis.

ESTIMATION OF TRACER "PLUME" CONCENTRATIONS Despite the absence of any measurements of the tracer-In concentrations in the flare "plume", we have a certain amount of information about its initial form and the mass of indium per unit distance along the burn track. Twelve flares were ignited simultaneously, six on each wing rack. Thus the rate of In-emission was 108.6 gm min-1 distributed equally between the two wing vortices (Dingle, 1968). The aircraft speed of 120 kt thus distributed 1.81 gm of tracer in each 0.0333 n. mi. (61.74 m.). The wing vortex may be envisioned as a cylinder of perhaps 2 m diameter initially. Thus the initial concentration X of tracer-In is 4.66 mg m3, and the total mass is 543 gm. Inasmuch as the tracer plume is a continuous cylinder, the dilution of the plume is entirely accomplished by lateral diffusion (longitudinal diffusion at the ends of the burn track is negligible). The actual dispersion, dividing, and transport of the material from the burn track is unknown, and is expected to be quite complex under the field conditions of 1 June 1970. To obtain an initial estimate of the airborne tracer-In concentrations prior to scavenging, we assume that the material was diluted by diffussion at the same basic rate as should prevail in dry air. The concentration, X, of a plume having source strength Q is (Turner, 1969): X(xrt) Q exp[-1/2(r/r) 2] (7) 2 7uc - 2 r where u is the mean axial (x-direction) wind speed, r is the 39

lateral dimension, and a is the lateral diffusion coefficient. For the present purpose, we assume that a =a a =a r y z under the convective situations we are dealing with, and that Q/u may be replaced by B/v, where B is the rate of emission of In from the burning flares and v is the airspeed of the airplane. Thus B/v is the quantity of In released per unit of distance along the burn track. The emission time of a sample of tracer-In, tb, is specified in terms of the distance, Sb, along the burn track, the air-speed, v, and the time of ignition, t: tb t + s /v. b o b Each element of the tracer plume at time t has therefore been diffusing for a time (t - tb), and the diffusion coeficient, r depends upon u(t - tb)* Values of -r have been published by Turner (1969), and by Islitzer and Slade(1968). There is uncertainty about the exact value of ar for our case, but the relative variation of a from one end to the other r of our 18-km-long streamer of tracer-In may be judged from Table 2. The figures are taken from Turner (1969). ho

TABLE 2: Values of a for beginning, ending, and mid-points r of tracer-In streamer at various times, t, (At = t - t u =.5 km min'l) tb = 1639 CDT 1644 CDT Average t CT (ignition) (burnout) ( midpoint) CDT uAt km a km uAt km a km uAt r 1644 2.5.475 1650 5.5.900 3.500 4.25.725 1700 10.5 1.60 8 1.30 9.25 1.45 1710 15.5 2.25 13 1.95 14.25 2.10 1720 20.5 2.85 18 2.55 19.25 2.70 1730 25.5 3.40 23 3.10 24.25 3.25 1740 30.5 3.95 28 3.65 29.25 3.80 1750 35.5 4.45 33 4.20 34.25 4.32 1800 40.5 4.95 38 4.70 39.25 4.82 Values of X along the axis of maxima are computed using the midpoint values of ar from Table 2. These are plotted in Figure 17 as estimated values of the air concentration along with the corresponding values of D. It is perhaps premature to attempt physical inferences from a comparison of X(x,o) against D, but certain features invite comment. Overall, X initially decreases in magnitude much more rapidly than D with increasing x(or t). There is, of course, a certain time interval between the airborne phase and the deposition, and this may account for some of the observed difference. The required ascent to accomplish nucleation, collection of droplets to form raindrops, and descent as rain to give the deposition appear to total to some 20 to 35 minutes of the average. The only process 41

of removal that does not require so large a time lag is that of fallout and washout of large particles. This process must exhaust the supply of large enough particles very rapidly, and so contribute only a negligibly small proportion of the deposition after a short time (10-15 minutes ~ 5 - 7.5 km). It is reflected by the early dip of D in Figure 16, but ^\ " x cannot be shown by the D values of Figure 17, nor should it be indicated in the X values because no particle-size separation is involved in the estimate of X( Aside from the above is the question of dispersion of the tracer from the wing vortices. We have as yet no information on the stability and persistence of these rotors, but we hope to get some help from airport studies done by Stanford Research Institute. In any event, the wing vortex behavior must prevent the rapid early dispersion indicated by the Sutton plume model to some extent. In addition the magnitude of a is uncertain as indicated above. Tetroon separation studies by Angell, Pack, and others (seeIslitzer and Slade, 1968) appear to give the closest approach to measurements of a at relevant altitude. The effect of the actual convective circulation field upon a remains one that probably cannot be measured directly. Finally, the imposition of the raining systems upon the field of dispersed material must have an effect upon that field which is not considered in deriving the X values, but which is present in the observed D values. It is not yet safe to say that all the difference between the rates of deposition and those of dilution is attributable to rain scavenging. Our modeling computations will eventually take account of the important processes and feedback L2

9 - 7- \ X(xo,t) (G/kM') 6543- \ 0 5 15 2025 0 5 10 15 20 25 X (KM) A X AND D ALONG STREAK AXIS FIGURE 17: Comparison of values of X and D along the deposition streak axis. 43

effects in the system so as to give a more quantitative view of the rain scavenging effect. LONGITUDINAL PROFILE OF DEPOSITION MASS A more comprehensive view of the deposition pattern is obtained by integrating the deposition across the pattern and plotting the result against downwind distance (Figure 2). The effect of the cross wind integration does not smooth the downwind variations in relation to those seen along the axis of maximum deposition values, D (Figure 16). This effect beyond 15 km appears to be related to the periodic growth, maturing, and decay of the convective systems observed in the radar and rainfall data. The pattern of variation of the integrated close-in deposition is quite similar to that shown by the study of the axial values out to about 15 km. The initial exponential decrease appears reasonable. The constant deposition region from 6 to 14.5 km is not expected, and requires further study. The data of Figure 2 provide a basis for the computation of the total deposited mass of In. The empirical curves and limits are indicated by dashed lines in the figure for which the equations are: (a) dM= 15 exp(-.105x) 0 = x < 6 km dM (b) = const 6 < x < 14.5 km (8) d (c) - = 17 exp(-.051x) 14.5 < x < 30 km -i~x."". "

Integrated within the indicated limits, these give respectively the following deposited mass totals: (a) M6 = 67 gm out to 6 km (b) M6_1 = 68 gm 6 km to 14.5 km 6-14.5 (c) M14530 = 91 gm 14.5 km to 30 km M = 226 gm total out to 30 km 30 which accounts for 42% of the total mass of indium generated by the flare burn within the sampling network. This excludes that part of the deposition collected in the extension of the pattern beyond 30 km, and it excludes all deposition south of the network. A residual mass of indium that is not at all shown by the deposition data is that which flows upward through the convective systems and out of their tops, and that deposited by processes that have not yet contributed within the 30 km limit (i.e., ice processes). The deposition between 6 and 14.5 km (Eq.8b) is envisioned as representing an overlapping effect of the processes represented by (8a) phasing out and the processes of (8c) phasing in. Thus the three domains of the deposition are (a) primarily fallout and washout of the largest particles; (b) the tail of (a) with the beginning of rainout of nucleating particles, the time delay being required for ascent of the nuclei and growth of the resulting cloud droplets prior to precipitation; (c) primarily nucleation scavenging with a diffusive collection component that increases its proportional contribution with downstream distance. The integration of (8) to infinity can be made to see how well the deposition accounts for the emitted tracer. 45

This integration accounts for 298 gm or 55% of the total of tracer-In. We can improve this accounting somewhat by estimating the deposition pattern south of the network. Subjective study of the patterns of Figure 1 in relation to the burn track suggests that the maximum of deposition is perhaps 2 km south of the network edge. The mass integral for this area, using (8a) extended to X = -2 km, gives 33.4 gm. The deposition still excluded from this region is that from the southern edge of the pattern in to the postulated maximum. Allowing a total of 45 gm deposition south of the network edge leaves a mass of 200 gm of indium for which no account is given by the present summation. INFERENCE OF DEPOSITION COEFFICIENTS The analogy with Engelmann's (1968) model is adopted here. In the present case, however, we do not have a widespread evenly concentrated tracer, so the concept of a reduction of concentration X over time t by scavenging processes in not applicable. Rather we propose to consider the total mass of tracer-In, I, as our basic point of reference, and consider the scavenging processes to act upon this quantity depleting it from the air and depositing it upon the ground. Hence It = I exp(-ct) (9) might express the total air-cleansing effect. Our data and our knowledge of the cleansing mechanisms, however, suggest that I is not a uniform population of b6

particles all having similar properties. Initially we assume that some fraction, fl, of I is in particles that will fall directly to earth upon release. We identify the following fractions: fl dry fallout fraction, r < 50 im, vt > 125cm/sec f2 impact-collection (washout), 2 < r < 50m f3 giant nuclei, 1 < r < 2pm f4 smaller active nuclei, 0.1 < r < lvm f5 "gap" particles, 0.03 < r < O.lim f6 brownian particles, r < 0.03pm Then dI =-I(E f +E f + f +..) (10) dt 1 1 22 33 expresses the total rate of removal, and the E. represent the respective partial scavenging coefficients. We assume the rate of removal equal to the rate of deposition, and transform the time dimension into distance by means of u, the mean speed of transport of the system. Thus - di dx -I( f + f +....) (11) dxr 1 1 2 2 I = I expl -(8ifl+62f2+-..)x u. and dI ~Elfl+2f2+-' (~ f +~ f +....)x I exp - f 1 (12) u u a) Initial Section: 0 = x < 6 km For the initial protion of the observed mass-deposition 47

profile (above), we have (8a): (IdM 1 =2 15 exp(-.105x) C x 1,2 The requirements for equivalence between (8a) and (12), assuming that dry fallout and washout are the dominant removal processes, are ~ f + f (a) 1 2 I = 15 gm km-1 u O f + f (b) 1 1 2 2 =.105 km-' u 15 Io = 143 gm o.105 and ~ f + E f =.0525 min1 1 2 2 This value for I is not satisfactory in view of the fact 0 that the emitted mass of In was 543 gm. (b) Second Section: 6 km < x < 14.5 km In this portion of the mass-deposition profile, dM dM _Z (8b) dx and we attribute this constancy to a compensation between the decreasing effect of fallout-washout and an increasing effect of nucleation. For such a compensation, the i. and the f. must be viewed as functions of distance (or its time 12I ~ ~ ~ ~ ~ b

equivalent). Thus the gross concept of these factors as constants, which is implicit in the integration leading to (12), is inadequate. We conclude that a satisfying evaluation can be made only upon a more detailed analysis of the processes that are glossed over by the above simplistic analysis. (c) Third Section: 14.5 < x < 30 km We have (8c): 17 exp(-.051x) IdM] dx 13,4,6 The requirements of equivalence between (8c) and X12) are E f + f + E f u f +Ef +Ef (b) 3 3 4 4 6 6 =.051 k (b) ~.051 kmu which leads to IO = 333 gm The approach of this value for I toward the known correct value, in contrast to the value deduced from the 0 = x < 6 km section of the data, suggests that the processes operating, as expressed by (e f + f f ), are much more nearly 3 3 4 4 6 6 constant in time than those operating in the initial section. Further study of the behavior of the Ei and f. is required to approach a complete understanding of the deposition pattern and processes. 49

REFERENCE Bleichrodt, J.F., J. Blok, R.H. Dekker, and G.J.H. Lock, 1959: the dependence of artificial radioactivity in rain on the rainfall rate. Tellus, 11, 404. Browning, K.A., 1964: Airflow and precipitation trajectories within severe local storms which travel to the right of the winds. J. Atm. Sci. 21, 634-639. Chem-Met Workshop, 1971: Verbal comment by a participant in the conference, held at Ft. Lauderdale, Florida, 25-29 January 1971. Dingle A.N., D.F. Gatz, and J.W. Winchester, 1969: A pilot experiment using indium as tracer in a convective storm. J. Applied Meteor. 8, 236-240. Dingle, A.N., and D.F. Gatz, 1971: Trace substances in Rain: concentration variations during convective rains, and their interpretations. Tellus, 23, 14-27. Dingle, A.N., and D.F. Gatz, 1971: Trace substances in Rain: Concentration variations during convective rains, and their interpretations. Tellus, 23, 14-27. Dingle, A.N., 1970: Rain Scavenging Studies, Report No. COO-1407-38, Contract No. AT(11-1)-1407, U.S. Atomic Energy Commission, the University of Michigan, Office of Research Administration, Ann Arbor. Georgii, H.W., and D. Wotzel, 1968: Spurenstoffkonzentrationem im Regenwasser in Abhangigkeit von der Tropfengrobe, Diplomarbeit. Newton, C.W., 1967: Severe convective storms. Advances in Geophys., 12, 257-308. New York, Academic press. 50

III. THEORETICAL MODELING 51

A BRIEF REVIEW OF ENTRAINMENT AND NUMERICAL EXPERIMENTS OF CUMULUS CLOUDS ENTRAINMENT MODELS FOR CUMULUS CLOUDS The first investigation of temperature and humidity distributions in and around trade-wind cumuli was the WymanWoodcock expedition to the Caribbean Sea (Bunker, et al., 1949). The in cloud temperature was found to depart slightly from that of the environment. Subsequent measurements of liquid water content have clearly shown that in most cases the liquid water contents in the main body of cumuli are usually only a fraction of those which would be expected from the assumption of themoist adiabatic accent of air. In order to overcome these shortcomings of the parcel theory, the concept of entrainment was introduced by Stommel (1947). There must be continual mixing between the rising parcel and the environment. This mixing of unsaturated environmental air into cloud is referred to as entrainment. The temperature difference between cloud and environment, and the differences between the observed liquid water content and the adiabatically predicated content, Wa, are the results of much mixing. i) The jet model Stommel (1947) noted that the temperature lapse rate in trade cumuli observed by Bunker et al. (1949) did not conform to expectations from moist-adiabatic. In order to explain this discrepency, it was assumed that at any level the ambient air is continuously drawn in through the sides of a cloud in a manner similar to that of a mechanically driven jet stream sucking environmental air into it. The drawing in of unsaturated air results in partial evaporation of the cloud liquid water content. This, combined with the 52

mixing-in of cool inflow from environment, causes the rising parcel to cool at greater than the moist-adiabatic rate. Stommel was inclined to the view that entrainment is due to turbulent lateral mixing, as in aerodynamic jet. His calculations of entrainment rate did not depend upon any such hypothsis. In this sense, his original treatment of the jet model was somewhat empirical. Austin and Fleisher (1948) gave a thermodynamic analysis of mixing and presented a general expression for the entrainment lapse rate. They demonstated that the lapse rate and the growth of a cumulus cloud are governed by the entrainment rate and the environmental dryness. They further suggested that such entrainment, called dynamic entrainment, is a necessary dynamic consequence of vertical acceleration. When vertical acceleration takes place, the vertical divergence must be balanced by a corresponding horizontal covergence. Unless the cloud contracts by at least a certain amount, entrainment is required to satisfy mass continuity. Houghton and Cramer (1951) hypothesized that the cloud is a steady-state and frictionless convective current, the diameter of cloud remains constant with height. Making the further assumption that the entrained air is uniformly mixed with the in-cloud air, they obtained equations for dynamic entrainment. It is clear that in Houghton and Cramer's treatment the entrainment rate is completely prescribed in advance by the assumption that the cloud cross-section does not change with height. As a result of evidence, the entrainment is size dependent. Haltiner (1959) presented a model similar to Houghton and Cramer's. In addition to systematic entrainment of 55

environmental air, he considered the transfer of heat and momentum by lateral turbulent diffusion. It is found that the effect of turbulent losses is the same as that of an exchange of matter between the cloud and its environment, the cloud losing some of its mass to the environment and acquiring an equal amount from it. In a study of trade-wind cumuli, Malkus (1954) considered that the momentum of the rising air is balanced by the momentum of ambient air entrained to satisfy temperature requirement and the momentum lost by detrainment. From continuity considerations, the net or dynamic entrainment is the difference between the gross entrainment and detrainment. The consistency of the various observations and computations gives some support to the picture of a cumulus cloud as a steady-state jet interacting with its environment by lateral entrainment and detrainment. Morton (1957) was the first to apply plume theory directly to convective clouds. He used equations for the conservation of mass, momentum, heat and moisture, and similarity considerations to deduce the radius, vertical velocity, temperature and moisture content of a plume entraining air from the: environment. Squires and Turner (1962) considered only entrainment and assumed its magnitude to be that suggested by Taylor (Morton, Taylor, and Turner, 1956) to compute the behavior of a thunderstorm updraft. They were able to compute the heights reached by clouds of known size in a given environment. This height was found to be a function of the assumed flux through the cloud base. Taylor's hypothesis states that for calculating entrain54

ment, the boundary of the jet (or bubble) may be considered as a surface over which the surrounding fluid flow in at a velocity proportional to the upward velocity at any level. The ratio of the entrainment velocity to the vertical velocity may be regarded as entrainment constant. Based on the results of laboratory experiments the value of a (entrainment constant) is about 0.1. According to Taylor's hypothesis the entrainment law can be expressed as 1 dM _2a M dZ = R where M denotes the mass flux, Z, the height, a the entrainment constant, and R the radius of the cloud. Thus the use _ dM of plume theory implies that the entrainment rate, ( = M dZ)' is inversely proportional to the cloud radius. As the horizontal dimensions of a cloud increases, therefore, the relative amount of environmental air entrained into it is reduced. In the foregoing studies, entrainment has been regarded as operative only from the sides of the cloud. Squires (1958) has suggested that dry air may penetrate above through the cloud apex. Probably so far no attempt has been made to provide a qualitative and quantitive estimate of which of these mixing mechanisms is more effective in their contributions to the distribution of liquid water content of cumuli. ii) The bubble model The entrainment idea of steady-state jet was challenged by Scorner and LudlUm (1953). A rising bubble or thermal consists of a ring vortex whose leading edges were eroded 55

away; the diluted air circulates around the outside and eventually is entrained into the wake. Scorer and Ronne (1956), Scorer (1957), and others used dimensional arguments and observational studies of fluid in stable stratified environments to show the detailed flow in and around a bubble. The bubble was observed to expand along a cone, its radius R being given by R = R + 3Z o where R is the value of R at the starting height (Z = 0) 0 and 3 is called the spreading coefficients. Various experiments indicate 3 f 0.2 to 0.25 irrespective of the buoyancy of the bubble. Levine (1959) has advanced a spherical vortex theory of bubble-like motions in cumulus clouds. He showed that the buoyant acceleration of convective bubbles is significantly diminished both by aerodynamic drag and by incorporation of still surroundings. On the basis of his picture of the bubble as a Hills Vortex, he showed that the magnitude of the entrainment can be written as 1 dM 9 K = M dZ 32 R where M is the cloud mass, R is the radius of the bubble, K is a dimensionless coefficient ranging from 0.5 to 0.8. It is reasonably well established that the entrainment rate of ambient air is inversely proportional to the updraft. Malkus and Williams (1963) demonstrated that the same form of vertical-rise-rate equation is applicable to both models. Various ways of measuring and reporting entrainment rate in convective clouds are compiled in Table 1 (Srivastava, 1964). 56

By applying jet theory to a large convective cloud with radius R = 5 km, and entrainment constant a = 0.1, the computed entrainment is 1 dM 2a 7cm M dZ - 4 x 10 cm. =M dZ R even for a large cloud, the entrainment cannot be neglected entirely. Indirect evidence suggested that strongest mixing is likely to be confined to a sheath surrounding the updraft core. There the ambient air which entrained into a large storm does not mix with the entire in-cloud air, only the outermost portion of the up-rising air is directly exposed to the eroding and dilulect effects surrounding flow. NUMERICAL EXPERIMENTS In the past decade, a number of attempts have been made to simulate convective clouds by nemerical experiments. Numerical models of convection in the atmosphere have progressed from the study of the motion of isolated dry thermals in an initially still surroundings to the development of clouds in vertically sheared air flow. Malkus and Witt (1959) were the first to apply numerical techniques to an isolated dry thermal model. It is assumed that motion is in a vertical plane with rectilinear symmetry. A lower neutral or unstable layer is topped by a more stable one, and a perturbation of potential temperature is introduced. They found that a mushroom-shaped bubble-like element was developed with the top part resembling a vortex-ring. Mason and Emig (1961) used a bubble model with two different formulations for rate of mixing with ambient air, both of which lead to the same conclusions. The important 57

parameters are the excess temperature and velocity of the bubble when it reaches the condensation level. Ogura and Phillips (1962) used a formal scale analysis to derive the equations of motion for shallow and deep convection in the atmosphere. The end result is a set of equations referred to as the "anelastic equations" If the vertical scale is small compared to the depth of an adiabatic atmosphere, the anelastic equations reduce to the Boussinesq equations. For shallow convection it is permissible to compute saturation vapor pressure without consideration of the effect of dynamic pressure on temperature, but the effects of changes in dynamic pressure on air temperature has a large influence on condensation rate in deep convection. Ogura (1963) developed a method to handle moist air (and later adopted by several other investigators). From the basic equations of hydrodynamics and thermodynamics, he derives the vorticity equation and conservation equations for total water substance and specific entropy of moist air. The model assumes a reversible process; that is, all condensed water is carried with the air and is available later for evaporation. Several computations with different lapse rates and ambient relative humidities are discussed. Orville (1965) used the numerical model of Ogura to study the development of cumulus clouds over mountains in an ambient wind. In the buoyancy term, virtual rather than actual temperature is used, and the consequences are discussed. Murray and Anderson (1965) assumes an Eulerian scheme for the field of motion, and a quasi-Lagrangian scheme for the fields of temperature, liquid water, and water vapor. 58

Evaluation of the fields of motion, temperature and liquid water is discussed in detail. In the models mentioned above, the precipitation, one of the important factors in the life cycle of cloud is neglected. Das (1964) presented a one-dimensional steady-initialstate model to study the effects of precipitation in convective clouds, the hydrodynamics and thermodynamics equations are solved with finite-difference techniques. He neglected the effects of entrainment of dry air into cloud as an initiation of the down-draft. He found that the start of precipitation initiates a sequence of events that results in the dissipation of the cell. Takeda (1966a) considered entrainment and precipitation. The cloud grows due to the release of heat by condensation until it reaches -20~C, at which time, raindrops are assumed to form. The falling of raindrops initiates the downdraft while the in-cloud temperature is still higher than that of the environment. Evaporation soon cools the downdraft. His work emphasized the downdraft near the cloud base. Takeda (1966a) incorporated the rainfall process into a two-dimensional model. Shear and precipitation are included. He assumes that initially a cloud has raindrops and horizontally uniform temperature on the upshear side, and only cloud drops and an excess temperature on the downshear Side. Thus a downdraft initiates - upshea; and an updraft - downshear. He derived a simple relationship between the rainwater content and the mean terminal velocity of raindrops. Srivastava (1967) presented a one-dimensional, nonentraining cloud model in which a simplified method of computing the mean terminal velocity of raindrops is used. Rain59

drops have exponential size distribution and grow by coalescence with cloud drops. Cloud water changes to rainwater at a prescribed rate when cloud water exceeds a critical value. Continuity equations for water substances include the production term which consists of autoconversions, accretion and evaporation (following the idea of Kessler, [1969] of parameterizing the cloud physics calculation). This type of treatment was extended to a two-space dimensional symmetrical model by Arnason et al. (1968). Liu and Orville (1969) use Srivastava and Kessler's techniques to study development of clouds over a mountain range. The treatment of evaporation of rainwater below the cloud is more general than Srivastava's. A parameter for the evaporation rate is introduced. Weinstein (1970) used a model very similar to that of Srivastava, but carried the work a bit further by considering entrainment, evaporation and freezing. The model is a onedimensional time dependent solution of the first-law of thermodynamic vertical equation of motion, and a series of moisture balance equations. The models mentioned above differ in many details, each emphasizing one or another aspect of the probelm, but nearly all of them have in common the property of one or two-dimensionality. Several existing models are based largely on Kessler's (1967) parameterized microphysics, the rate of conversion from cloud water content to rain water current was assumed and cloud physical processes were parameterized. 60

REFERENCES Arnason, G., R.S. Greenfield and E.A. Newburg, 1968: A numerical experiment in dry and moist convection including the rain stage. J. Atmos. Sci., 25, 404-415. Austin, J.H., and A. Fleisher, 1948: A thermodynamic analysis of cumulus convection. J. Meteor., 5, 240-243. Bunker, A.F., B. Haurwitz, J.S. Malkus and H. Stommel, 1949: Vertical distribution of temperature and humidity over the Caribbean Sea. Papers in Phys. Ocean. and Met., MIT and Woods Hole Ocean. Inst., 11, No. 1. Das, P., 1964: Role of condensed water in the life cycle of a convective cloud. J. Atmos. Sci., 21, 404-418. Haltiner, G.T., 1959: On the theory of convective currents. Tellus, 11, 4-15. Houghton, H.B., and H.E. Cramer, 1951: A theory of entrainment in convective currents. J. Meteor., 8, 95-102. Kessler, E., 1969: On the distribution and continuity of water substance in atmospheric circulations. Meteor. Monogr., 10, No. 32, 84 pp. Levine, T., 1959: Spherical vortex theory of bubble-like motion in cumulus clouds. J. Meteor., 16, 653-662. Liu, J.Y., and H.D. Orville, 1969: Numerical modeling of precipitation and cloud shadow effects on mountain-induced cumuli. J. Atmos. Sci., 26, 1285-1298. Malkus, J.S., 1954: Some results of trade-cumulus cloud investigation. J. Meteor., 11, 220-237. Malkus, J.S., and R.T. Williams, 1963: On the interaction between severe storms and large cumulus clouds. Meteor. Monogr. 5, No. 27, AMS, Boston, Mass., 59-64. Malkus, J.S., and G. Witt, 1959: The evolution of a convective element. A numerical calculation. The Atmosphere and the Sea in Motion, New York, Rockfeller Inst. Press, 425- 39. 61_

Mason, B.J., and R. Emig: Calculations of the ascent of a saturated buoyant parcel with mixing. Quart. J. Roy. Meteor. Soc., 87, 1961, 212-222. Morton, B.R., 1957: Buoyant plumes in a moist atmosphere. J. Fluid Mech., 2, 127-144. Morton, B.R., G.I. Taylor, and J.S. Turner, 1956: Turbulent gravitational convection from maintained and instantaneous sources. Proc. Roy. Soc. A., 234, 1-23. Murray, F.W., and C.E. Anderson, 1965: Numerical simulation of the evolution towers. Report SM-49230, Douglas Aircraft Company, Inc., Santa Monica, Calif., 97 pp. Ogura, Y., 1963: The evolution of moist convection elements in a shallow, conditionally unstable atmosphere. A numerical calculation. J. Atmos. Sci., 20, 407-424. Ogura, Y. and N.A. Phillips, 1962: Scale analysis of deep and shallow convection in the atmosphere. J. Atmos. Sci., 19, 173-179. Orville, H.D., 1965: A numerical study of the initiation of cumulus clouds over the mountainous terrain. J. Atmos. Sci., 22, 684-699. Scorer, R.S., 1957: Experiments on convection of isolated masses of buoyant fluid. J. Fluid Mech., 2, 583-592. Scorer, R.S., and C. Ronne, 1956: Expdriments with convection bubbles. Weather, 11, 151-154. Scorer, R.S., and F.H. Lundlam, 1953: The bubble theory of penetrative convection. Quart. J. Roy. Met. Soc., 79, 94-103. Squires, P. 1958: Penetrative downdrafts in Cumuli. Tellus, 10, 381-389. Squires, P., and J.S. Turner, 1962: An entraining jet model for cumulonimbus updrafts. Tellus, 14, 422-434. Srivastava, R.C., 1964: A model of convection with entrainment and precipitation. Scientific Report MW-38, Stormy Weather Group, McGill University. 62

Srivastava, R.C., 1967: A study of the effects of precipitation on cumulus dynamics, J. Atmos. Sci., 24, 36-45. Stommel, H., 1947: Entrainment of air into a cumulus cloud. J. Meteor., 4, 91-94. Takeda, T., 1966a: The downdraft in the convection cloud and raindrops: A numerical computation. J. Meteor. Sci., Japan, 44, 1-11. Takeda, T., 1966b: Effect of the prevailing wind with vertical shear on the convective cloud accompanied with heavy rainfall. J. Met. Soc., Japan, 44, 129-144. Weinstein, A.I., 1970: A numerical model of cumulus dynamics and microphysics. J. Atmos. Sci., 27, 246-255. 63

ON THE VERTICAL DISTRIBUTION OF CONTAMINANT - THE MODEL INTRODUCTION It is well known that the main mechanism of removing contaminants from the atmosphere is their rainout and washout. Greenfield (1957) proposed a theory of rainout and washout of aerosol contaminants from the atmosphere. A munber of models of the same general type have been proposed. Much attention has been given to the processes of removing contaminants, but the input mechanism and the initial vertical distribution of contaminants in a severe storm have been investigated very little. Gatz (1966) calculated the input rates of radioactivity and of moisture to the storm using a procedure adapted from Newton and Fankhauser (1964). It is understandable that the problem of the distribution of a contaminant in a big cloud is complicated by dynamic features. The primary interest of the present investigation is to construct a tentative model to explain the initial distribution of contaminants and the input mechanism. THE MODEL The updraft is modeled as a circular tube having nonuniform properties both in radial and vertical directions. Entrainment is taken into account. The horizontal cross section radii are allowed to vary with the height of the storm. The contaminants transported into the storm are considered mainly due to the inflow of substance up through the cloud base. The cloud is assumed to be a steady state 64

system moving with a horizontal speed. THE UPDRAFT RADIUS In the present study, a cloud core is assumed to be the area occupying an updraft with a profile specified by the following analytical expression W = W(r,Z) = 0.91 Wmax cos 2R(Z sin 7 - sin H 1) H 4 H where R (a function of height only) and H are the radial and the vertical boundaries of the modeled storm, respectively, and W(r,Z) is the vertical velocity. With this profile, we will have the following characteristics of updrafts in cloud a) the updraft increasing with height having a peak at the level Z = 0.62 H (Fujita and Grandoso, 1968), b) a gradual decrease in W with distance from draft center to the draft periphery. The form of the horizontal variation is considered to be axially symmetric so as to maintain mathematic simplicity of the problem. To consider the eroding effects on the velocity profile at the edge of a draft by mixing due to entrainment, the horizontal profile is assumed to be wave shaped. The mass flux through a horizontal circular cross section perpendicular to W can be written as R(Z) M (Z) = p(Z) 27rr W(r,Z) dr (2) 10 where P (Z), the air density, is assumed to be a function of h eight only. 6.

Substituting Eq. (1) into Eq. (2), we lave TrZ 1 2TFZ R(Z)'irr M(Z) = 0.91 W 21p(Z) sin s in os 2 rdr. max R 4 H 2R(Z) This can be written as M(Z) = 3.64 1 - - W R2(Z)p(Z) sisi n sin (3) max 4 H' We differentiate Eq. (3) with respect to Z after taking the natural logarithm and introduce the function rZ 1 2TrZ *~ W(Z) = sin sin. - 4 H Hence 1 d (Z) 1 d R2 (Z) dp(Z) 1 dW(Z) M (Z R W(Z) + + (4) M(Z) dZ R(Z) T dZ p(Z) dZ W(Z) dZ It is reasonably well established that the rate at which environmental air is entrained through a vertical cloud wall is inversely proportional to the draft diameter and can be written as 1 dM(Z) _ 2a M(Z) dZ R(Z) (5) where H is the entrainment rate and a is entrainment constant, determined in tank experiments a = 0.1. Combination of Eq. (4) and Eq. (5) gives dR(Z) R(Z) d dZ + 2 [ln p(Z)W(Z)] =a (6) with the initial conditions at the cloud base Z = ZO, p = Po' R = Ro, and W = W. The solution of Eq. (6) is expressed 66

as =^"^ ^T W[p()WZ) " Ep(Z)W(Z)] dZ (7) R(Z) R (' Z)]P o 1/2 + c(Z)W(Z) l/2f [p(Z)W(Z)] d 0 THE UPDRAFT TRAJECTORY Bates (1961) devised a numerical method to compute the trajectory of an updraft. A draft is treated as if it were composed of a stack of wafers with uniform draft velocity in each wafer, and entrainment was neglected. A horizontal wafer is subjected to an aerodynamic drag force resulting from the relative motion of the ambient wind field and the wafer. Bates and Newton (1965) computed the two-dimensional trajectories of wafers rising from the cloud base by assuming a vertical velocity of the wafer and a constant horizontal storm velocity. Thus, they found that the shapes of the updraft profiles are very sensitive to the value of storm velocity. To apply the principle (Bates, 1961; Bates and Newton, 1965) of the computation to our model, we consider that the mass of the wafer varies due to the fact of entrainment. The drag force in a slab of thickness AZ of the cloud can be stated as F= Ce (Z)V R(Z) AZ D D e r where CD, the drag coefficient, is taken to be 1.8, and V = D' r (V - V ) where V(Z) is the ambient wind speed and V (Z) the w w horizontal speed within the wafer. This may be equated to 67

the rate of the change of momentum d d (MWw) = CDPe(Z) V R(Z) AZVr where M(Z) = TR2(Z)p(Z)AZ is the mass of the wafer. In a steady state, 1 dM - dM _W dM W M 2a M dt M dZ M Z W R(Z); similarly dV dV dV w w dZ w dt dZ dt dZ Thus dV (Z) dZ V 2~"z3) ( 8 ) dwZ =a- a 2Vw(Z) + a3V2(Z) (8) where where ~CDV2(Z) a1 =7(Z) 2V(Z)CD 2a D 2 20 2 =TR(Z)w CD a = 3 7TR(Z)w In our cloud model, the updraft is non-uniform across the circular cross-section, to simplify the computations, an area mean is defined as 27; R(Z) A = 2 R() Ar dr R2 (Z) J 68.

Hence the area mean updraft WH at any level may be written as [using Eq.(1)]. 3.64(7T-2) W Z sin 1 27 W - l~ ia wx { sin sin ] H 72 max H 4 H The environmental wind field V(Z) is based on Bates and Newton's (1965) profile, varying linearly from -5 m/sec at the cloud base to 30 m/sec at 8 km and decreasing above that. V(Z) can be written as 35 V(Z) = -5 + 8000 I Z Z > 8000 m and V(Z) = 30 - ( Z 8000) Z > 8000 m Eq.(8) may be integrated analytically by assuming mean values of all parameters for a thin layer. Considering first the case where A = 4a a - a2; A > 0 1 1 3 2 1 we can then write the solution of Eq. (8) in the form ( 2aV - a /Av Vj+l [a + tAn tan tanI 3w + 2- (Z.j-Z) ]/2a (9) w 2 1 2 j+J j 3 1 Turning now to the case where A < 0, A = a2 -4a a 1 1 2 1 3 69

We write the solution to Eq. (8) in the form V j = ((C + a) + (/r - a)2 - expA (10) where 2a V - a- A 3 w 2 A2 = A (Z Z) + log 3j -2 a j+l 2aVIaV]/a 3 w 1 the time required for a wafer to ascend from the layer j to the layer j+l is expressed by At= i dt= IZj+ (Z Z o Zj. W W 0 J Where W is the mean vertical velocity of the layer. From a knowledge of the horizontal speed of the draft at all levels, the horizontal displacement in a relative frame fixed to the cloud moving with speed C can be computed from N N X(Z) = Xj(Z) = [V (Z) - C](Zj+ - Zj)/ j=1 j=l w THE CONTAMINANT DISTRIBUTION Consider a thin slab of thickness AZ in a cloud, the total content of a tracer element, say Q, is Q(Z) = TR2(Z) p(Z) q(Z) AZ (11) where q(Z) is the mixing ratio, a function of height only. 70

Since all terms on the right hand side of Eq.(11) are functions of height only, if a steady state is assumed, then we differentiate Eq. (11) with respect to t (time) 2p R21 NW ( dt V= TrAZ7qW q( + Rq -z ++ pR2q -z + P2 (12) when entrainment is taken into consideration, for a balanced state, the total rate of change of tracer element Q must be equal to the lateral entrainment, which can be written as [dQ1 27R(Z) u P (Z) q (Z) AZ where u is the lateral entrainment velocity, P (Z) is ambient e air density, for simplicity, P(Z) = P (Z), q (Z), is the tracer element in the environmental air. For balance C dQ dividing by R2(Z)Wq(Z)P(Z) and assume = 0.1 (Squires and Turner, 1962) we have W ( p 1 aR2) + 1 (Z) 0.2q (Z) = R(z)q(z) q1 R2(Z) But 1 3p(Z) 1 a2(Z) 1 W 2o p4Z) z (Z) + = RT(Z' [Eqs.(4) and (5)] R2(Z) zaz az 2aq (Z) q(Z) 2 q(Z) 2 (14) ~-z + R(Z) R(Z) For reasons of simplicity, R(Z) determined from Eq.(7) may be considered as a constant at any level of interest. 71

The solution of Eq. (14) is obtained by putting the initial conditions at the cloud base Z = 0, q(Z) = q (0). Thus, we have 2a 2 O 2z q(Z) = q(O) e R(Z) 2a R(Z- R(Zq(Z) = q(0) e Re+ R( ) qe(Z) e dz (15) If we assume that the total amount of tracer element Q in a slab is non-uniformly distributed and the horizontal distribution function of mixing ratio q has a wave form q(r,Z) = A cos 2R() q(Z) (16) where A is the amplitude of the distribution function. Then Q(r,Z) = A 27 P(Z) q(Z) cos 2R( r dr AZ (17) L'0 2R(Z) Integration of Eq. (17) gives Q(r,Z) = A 4(2) [P(Z) R () q(Z) AZ] (18) From Eq.(11) and Eq. (18), we obtain 2 7T2 A = 4(-2) 74 (-2) and q(r,Z) = 4(2) q(Z) cos 2R(Z) (19) 72

Hence the final expression for mixing ratio is found by combining Eqs. (15) and (19) iT2'r z2 Tr- 2a R-Z) q(r,Z) = 4'2) cos 2RZ) q(0) e R(Z) + (Z) e 4^ ^ 40T-2) ^ 2R(Z) RLUZ e' (2)e (20) z (Z) e() dZ' 0 COMPUTATIONS AND RESULTS For the purpose of obtaining the features of contaminant distribution within a large severe thunderstorm, we selected certain parameters which remain unchanged throughout the computations. The cloud domain under consideration has a depth of 10 km. The entire cloud depth was divided into layers of equal thickness, AZ = 100 m, extending from the base at Z = 1.5 km to the top of the cloud at Z = 11.5 km, thus establishing 100 layers between the cloud base and the top. A constant lapse rate of 0.650C/100 (an unrealistic assumption, but does not nullify the qualitative conclusions) with a surface temperature T = 2880K was used to obtain both in-cloud and ambient air density functions. Figure 1 shows the vertical velocity profile at the center of the cloud with a maximum velocity Wa = 35 m/sec, max the maximum velcoity occurs at H = 7.1 Km above ground. The variation of mean draft radii with height is illustrated in Figure 2. In high levels there is considerable broadening due to the rapid decrease in wind velocities. At 9.5 Km (above cloud base), the draft radius R = 10.6 Km 73

1 9 8 7 6, 5/ 45 bo 4_ 3 2 1 5 10 15 20 25 30 35 40 Vertical velocity (m/sec) FIGURE 1: Verticle wind profile. 74

10 9 L 8 6 5 _ 4 43 3 2 1 1 2 3 4 5 6 7 8 10 11 12 13 1 Mean Draft Radius (Kh ) FIGURE 2: Mean draft radius. 75

is approximately twice as large as the initial cloud base radius R = 5.0 Km. Minimum draft radius R = 3.25 Km occurs at Z = 4.1 Km (above cloud base) which is about 1.5 Km below the level of the maximum wind. The horizontal velocity in a draft with varying radii illustrated in Figure 3 is computed from Eq. (9) for V(Z) < 0, and from Eq.(10) for V(Z) > 0. The draft profile corresponding to the case Wmax 35 m/sec, C = 27 kt = 13.9 m/sec is shown max in Figure 4. The horizontal displacement of the draft element decreases upward with a maximum value of about 573 m in the first layer and a minimum value of about 18 m at the top layer, the total horizontal displacement for a wafer rising from the cloud base to the cloud top is approximately 12.7 Km, in this ascent, total time required is about 18 minutes. As would be expected, the stronger the vertical velocity, the smaller is the time required for the rising air to travel through that layer, and the smaller is the horizontal displacement during the time of the wafer's ascent. As pointed out by Bates and Newton (1965), the shape of the updraft varies considerably according to the translational velocity of the cloud system. In our computations, the storm speed C = 13.9 m/sec > V (max) = 13.4 m/sec has been assumed to be that which makes the updraft lean into the wind, which is generally considered as a realistic condition through most of the cloud depth. Figures 5 and 6 show draft configuration based on different storm speeds. Mixing ratio in the vertical dimension is evaluated from Eq. (15) under the assumption that the environmental mixing ratio of contaminant is zero, and the initial value 76

10 - C=8 m/sec - 9 8 - ^ T storm speed 7C 13.9 m/sec _ 7 %MOO^~~~~~~~~~~~~~ ^ ~(ENV.) w (wafer speed) Sri/ 4 A ^IL Wmax = 35 m/sec 4,m 3 2 F / | / l l SAC = 20 m/sec -5 0 5 10 15 20 25 30 (m/sec) FIGURE 3: Horizontal wind speed.

10 F^. ^-z ~draft axis - draft / s8 -radius J ^~\ \ 4 draft boundary?jb~~b 7 3r storm speed C = 13.9 r/sec 2~~~~~~~~~~ 25 20 15 10 5 0 3 Horizontal Displacement (1(m) FIGURE 4: Draft Shape.

10,^draft axits 9~~~ draft 8 radius 7 6~ / d^~ I^~~ ^ ^2~?-draft boundary 5 bO 4 storm speed o i C 8 cm/sec 3 2 " 1^ \ 20 15 10 5 0 5 Horizontal Displaciment (Km) FIGURE 5: Draft Shape

10 10^ N^^-~draft axis 9 draft / radius 7- \ \cc I 7 6g~ \A-0'v0%.~~~~2 draft boundary 45 4) storm speed 3 _C 20 m/sec 2 1 30 25 ~ 20 15 10~ 5 Horizontal Displacement (Km) FIGURE 6: Draft Shape.

at cloud base is set at unity. In the case of R = 5.0 Km, the mixing ratio decreases from 1.0 at Z = 0, (cloud base) to 0.62 at Z = 9.5 Km. In practice, the contaminant content Qi is computed and summed over discrete layers from the cloud base to the height 9.5 Km using the relation N N Q = Qi = 1 TR (Z) p (Z)q (Z)Az, N = 95 AZ = 100 m i i=l i= where i is the layer index and N is the number of layers in the cloud. The presence of a bar over a parameter indicates thet its mean value in the layer is used. qi(Z) is determined from Eq. (15). Q is the total amount of tracer element in the cloud (if we assume that no air is detrained from the cloud). One obtains the relative amount of contaminant in each layer from the realtionship Q Qrel^ (21) rel = The profile of Eq.(21) is shown in Figure 7. It is found that the minimum content occurs at the maximum vertical velocity region in the layer from 5.6 to 5.7 Km. The mean mixing ratios of each layer in the vertical dimension are converted into units of gm of contaminant per gm of the air from the relationship q(Z) = rel f amount of contaminant g T R(Z) (Z)AZ amount of air g m/gm ( the variation of qi(Z) with height is shown in Figure 8. 81

10 0-oo 9 90 8 -80 7 0to 6 -60= 05 — 4 -40 ij 1r 1\ 0 10 -3. 1.0 2.0 3.0 Mass of Contaminant (% total mass of contaminant) FIGURE 7: Vertical profile of contaminant in a disc of thickness AZ = 100 m. 82

10 100 9. 90 8 80 7 _70 6 60 S 5-z 50 -P F \ 3 fi \ 4 40 3 302 20 1 - 10 2.0 3.0 4.0 Mean Mixing Ratio (10-15gm/gm) of Tracer FIGURE 8: Vertical profile of mean mixing ratio in a disc of thickness AZ = 100 m. 8~

10 9 8?7 (1) 6 2) X (1) Base Layer 1.5-1.6 Km f% j ) \ ((2) Mid Layer 6.5-6.6 Km k 5 _ H\\ (3) Top Layer 10.5-10.6 Km 3. 5 1 2 3 4 5 6 8 9 Mean Draft Radius (Km) FIGURE 9: Radial profiles of mean mixing ratio. 84

Mass of Contaminants (% of total contaminant mass) HCI eni o::y H o (D 0 o0 0 0 0 0 ^M ~~I ~d ~u ~-'~~ ~ ~~ I\) pTT ~J Q1 - ^ n~ 1 \ FI.I___-i,. \ MD ~U!"^ "I *~~ ~, p.l ^ "~ CD a~~ ~ ~ ~ ~~~~~~~~~ st ( ^0^ r o ^ aI - Qpi 1 -t, 0- - "^ u ^ -I ~ LCn Ii C: opi H la, 0 0 I-h ~ ~ ~ -~ (D~~~~~~~~~~~~~~~~~~~~~~~~D C cm x0 OM 0 ^~D0 COO H-* F-i ^ \ ~ i <~~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~ ~~~~~~~H U- ^^\ r M fl) H~- 0 @ HD ^X 2 ~~~~~~~~~~~~0 (D H U) FJ -J0 (n p-~ ui * ^ ~ ~~~~~~~~~~~~00^ D"~~~~~~~~~~~~~~~~~~~~~~~~( ^' ~ ~ ~~~~~~~~~~

The horizontal distributions of mixing ratio in certain layers are computed from Eq. (20) by substitution of qi(Z) from Eq.(22) and illustrated in Figure 9 which shows the amplitude of mixing ratio decreasing with height, and considerable narrowing in middle levels and broadening in high levels. The mass of contaminant (percentage of total mass of contaminant in the cloud) in each cylindrical shell of thickness 500 m (extending from cloud core to the cloud edge) and of thickness Z = 100 m is computed from Eq. (17) using mean values of all parameters. The distributions of mass illustrated in Figure 10 present a wave form with peaks approximately in the middle shell. CONCLUSION As mentioned in the introduction relatively little is known about the actual tracer distribution within a larger thunderstorm. We offer the preliminary results of the oversimplified model of numerical experiments aimed toward providing a general picture of the distribution. Our intention is to examine different processes (Nucleation, etc.) separately and to introduce them into the model. We recognize that a model of this nature is, of course, far from conclusive. But it is hoped at least that tracer estimates obtained in this simple manner will be found useful, until more and better measurements have become available. 86

REFERENCES Bates, F.C., 1961: The great plains thunderstorm - a model. Ph.D. dissertation, St. Louis University, 164 pp. Bates, F.C. and C.W. Newton, 1965: The forms of updrafts and downdrafts in cumulonimbus in a sheard environment. Paper presented at American Meteorological Society Meeting, Reno, Nev. Fujita, T. and H. Grandoso, 1968: Split of a thunderstorm into anticyclonic and cyclonic storms and their motion as determined from numerical model experiments. J. Atmos. Sci., 25, 416-439. Gatz, D.F., 1966: Deposition of atmospheric particulate matter by convective storms: The role of the convective updraft as an input mechanism. Ph.D. dissertation, University of Michigan, 216 pp. Greenfield, S.M., 1957: Rain scavenging of radioactive particulate matter from the atmosphere. J. Meteor., 14, 115-125. Newton, C.W., and J.C. Fankhauser, 1964: On the movements of convective storms, with empahsis on size descrimination in relation to water-budget requirements. J. Appl. Meteor., 3, 651-668. Squires, P., and J.S. Turner, 1962: An entraining jet model for cumulonimbus updrafts. Tellus, 14, 422-434. 87

CLOUD AND RAIN FORMATION INTRODUCTION The purpose of this study is to apply a numerical model of convective cloud, somewhat different from the models presented previously, to the scavenging process by rain. In this study, we attempt to relate the cloud model to the scavenging model. The first part includes the cloud model, and the second part includes various scavenging processes. We restrict the problem to one space dimension, the veritcal, and to time. Any model that can provide a correct picture of the main vertical currents should not be very far from true. This is the view taken in the present investigation. If the model can somewhat explain the experimental data then it can be extended to a more elaborate model. ASSUMPTIONS The following assumptions are used in the model: 1) The model is time dependent, one dimensional in the Z plane. 2) The terminal velocity of water drops (Gunn and Kinzer, 1949) is described by using a regression equation (Dingle 1970). This equation is used through this study. 3) The size distribution of water drops follow an exponential function of diameter (Marshall and Palmer, 1948). 4) It is assumed that the cloud remains just saturated with respect to liquid water, evaporation will be considered to occur only from cloud drops in the entrainment process, in the subcloud region, raindrops may evaporate. 88

EQUATION OF MOTION The equation of motion used in the model may be written as d_(MW) = (TwJ gM or and Tv = (1 + 0.61 q5)T. T (1 + 0.61 qe)T. The buoyancy force per unit mass is[(T - T )/T - Qw]g v ve v w where (Tv - Te) is the excess virtual temperature of the parcel over that of its environemnt at the same level. The liquid water mixing ratio Qw(gm/gm), which makes a significant contribution to the buoyancy, is included in the equation. THE ENTRAINMENT The entrainment lu in tis included in this model cloud because it cannot be neglected entirely even for a large cumulonimbus cloud. The theories imply that the entrainment rate v is inversely proportional to the cloud radius and can be expressed as 1 dM 2a M dz Rc 89

where M is the cloud mass, z is the vertical distance, Rc is the cloud radius and a is the proportionality constant, from laboratory experiments on plumes it is found to be close to 0.1. dM dM dt Since Sine dz dt dz 1 dM thus W = dt EQUATION FOR WATER SUBSTANCE The liquid water mixing ratio, Qw' of the modeled cloud is assumed to be the sum of cloud water mixing ratio Qc' and rainwater mixing ratio Qr. Thus Q = Q + Qr (gm/gm). The cloud water mixing ratio is defined as the water substance of very small droplets which have a negligible terminal velocity and thus are carried along by the veritcal current. Due to continuous condensation, collision and coalescence amongst themselves, some of the cloud droplets grow large enough to have an appreciable terminal velocity, which can fall against the updraft as rainfall. The water substance of these larger drops is defined as rain water mixing ratio. The rate of change of water substance of a rising air parcel in which there is entrainment is assumed to be balanced by the water substance in the entrained air. The equation reads d dM dt [M (q + Qw) ] d- qe 90

Expanding the above equation leads to dQw 1 dM dt =dt - e wM dM Following Kessler (1969), the rainwater mixing ratio changes by three processes representing a) the accretion of cloud water by rain drops, b) the conversion of cloud droplets to raindrops, and c) the evaporation of raindrops in the subcloud regions. a) Change of rainwater mixing ratio through accretion of cloud droplets. The rate at which volume is swept out by a single raindrop of diameter D falling at speed Vt is iD2/4 x Vt, and the raindrop coalesces with the cloud droplets in its path with an efficiency E, thus the growth rate of a single drop is dm Q TrrD2VtE dm) = t dt a 4 The total rate of increase of rainwater mixing ratio by the capture of cloud drops is given by dQ Q TD2V E (r-) = N dD. dt a 4 D If it is assumed that rain drops of terminal velocity equivalent to or in excess of the updraft will fall to collect the cloud droplets, and the raindrops of diameter greater than a certain value will break up, then the integration limits of the above equation is from D = D(Vt > W), D = Db (breaking up diameter) rather than from D = 0 to D = o. Then dQr Db Qc TDTD2E dt ),a.....c.. NDdD. -.-.(3) dt a 4 D c,/ 91

Evaluation of the integrals occurring in the above equation, Vt, ND and E must be specified. Gunn and Kinzer (1949) made accurate measurements of theterminal velocities of water drops. Their observations were made at a pressure of 1013 mb, a temperature of 20~C, and a humidity of 50%. In order to apply the experimental data to the drops at different altitudes, it is necessary to make some reasonable estimates. Using statistical techniques, Dingle (1970) found an empirical relationship between terminal velocity and drop's diameter, which can be written as V = A + A D + A D2 + A D to 0 1 2 3 with Ao = -16.603 cm/sec Al = 491.844 sec-1 A2 = -88.801 cm-lsec- A3 = 5.488 cm-2sec-1 where V in cm sec-1, D in mm to an approximate correction is applied to the empirical equation to evaluate the terminal fall speed at different altitudes. This correction can be written as Vt z Po 1/2 t,o where V is the terminal spped at air density p and V t,o 0 t,z that at air density p. Observations by Marshall and Palmer (1948) of rain at earth's surface have shown that the drops, on the average are size-distributed in approximate accord with an inverse exponential law. We assume that the raindrop size distrihas the same form throughout the cloud. 92

-XD ND = Noe where No and X are parameters of the distribution. Introducing these expressions into Eq. (3) we obtain dQ, NOf pQ P, b - (dQr) NOQC( )1/2 (AoD2 + AD3 + A2D + A3D )e EdD. dt a 4 P D c Integrating the equation, if E is independent of D, we obtain m dQr No07TQcE 1/2 i -XDb -r ^r, Po 1/2 f (_o)r' -XDb m-r -)' A e (e1 D... dt 4 P n rO (mr)!(X) b e DcDm-r)]-(4) where m = 2, 3, 4, 5 and n = 0, 1, 2, 3, and m = n + 2. The total rainwater mixing ratio, Qr' is given as T b -XD Q I NoD'e' dD. Qr = 6 D c As before, solving the integral gives 6Q X4( 0) + 3(De-XDb - D3e-DC) + 3X2(Dbe-Db D2eDc) p b c b c + 6X(Dbe -Db DeXDc) + 6(e-XDb -e c) - - be D C) = e which given the relationship between X and Qr The cutoff diameter. D can be calculated with the help c of the empirical equation V =A A + AD) + A3D =c W (updraft) t,z (A~ + ADc D C 3C) ( po) 1 /29 93

and Db can be determined from the experimental data. The air density p is assumed to decrease with height as e (g/RYe 1) P Po T eo b) Cloud Conversion According to Kessler's argument, a threshold value, Qth for cloud water mixing ratio is specified. As Qc' cloud water mixing ratio passes this threshold value, the rate of cloud conversion increases with the-cloud water mixing ratio, but depends on the conversion factor a.:Below the threshold value, no cloud water is converted. Such a process is defined mathematically as dQr { a(Qc - th) Qc th "dt c (5) ~ Qc <th c) Evaporation of Raindrops The evaporation for a single drop is dm 2 TD(S-l) (1 + FR1/2) (gm/sec) (dt^e =-A+ ^Be (gm/sec/ ) (Byers, 1965, dt e A + B e p. 113) where S is q/q, the ventilation factor F = 2 for R > 5.6 s e (Abraham 1968), and Re, the Reynolds number a function of drop sizes. A and B are a latent heat term and a diffusion term respectively. The rate of change of rainwater mixing ratio due to 94

evaporation is given by dQr) = (2(S) N D(1 + FR/2)e dD (6) dte p(A + B) D e c where we assume (S-1)/(A + B) is drop size independent, but R is a function of drop diameter, and e VtD R = _ e "k Introducing these expressions into Eq. (6) we obtain dQ 2 I2T(S-1)N dQr) = 0~rS1N 1 (De-DC D e-XDb) + 1/X2(e.D dt e p(A + B) (D - b ) + (7) e0 1 2 -D D e D e XDb) + FA +2 A3 D3 eD D zVk C If it is assumed that,in the cloud, only cloud drops will be considered to evaporate in the entrainment process, below the cloud base, raindrops may evaporate depending on the vapor mixing ratio. With the help of Eq.(2) which gives the relationship for mixing ratio, and the various conversion terms given above [(4), (5) and (7)3, the equation for water substance can be sumarized as dQw dQ dQr dt dt dt d dQs 1 dM 1 dM dQc d c dt —" -M d te Md WFc+() (") - c(- ) (8) 95

On the right side of Eq. (8) the first term represents the production of Qc by condensation, the second term represents the evaporation of Q in the entrainment process, the third term represents the dilution of Qc, the fourth term represents the accretion of Qc by Q discribed by Eq. (4), the fifth term represents the conversion of cloud droplets to raindrops. Similarly, the continuity equation of Qr is dQ 1rdM dQ dQr dQr 1 dMrr dt M d r dt a + dtc dt e' where the first term represents the dilution of Q by entrainment of clear air and dq _ dq dQr dtdt t dt e dQ dQ dQ dQc r, ^c r c dta dt ) a' dt ) c dt THERMODYNAMICS We consider that the parcels include dry air, water vapor, and condensed particles. The cloud air is assumed to remain just saturated. In the convective current considered here, there is entrainment of unsaturated air which becomes saturated by evaporation of some of the condensed moisture. The equation for cloud temperature is very similar to those outlined by Austin and Fleisher (1948), Haltiner (1959) and Weinstein (1970). As the parcel ascends, a certain amount of water vapor, dqs, condenses. The heat released through condensation is dQc = -LMdqs. 96

While rising, the parcel entrains non-saturated environmental air. The originally saturated cloud air becomes slightly subsaturated, in order to saturate the entrained air, heat is required to evaporate some of the liquid water. This heat loss is a latent heat transfer and is given by dQ = -L(q - q)dM. The heat needed to raise the temperature of the environmental dry air from T (its original temperature) e to T (the cloud temperature) is dQ = -C (T - T )dM. s p e Whenever raindrops fall into an unsaturated region, for example, below the cloud base, the water is evaporated until the wet bulb temperature is reached. The heat loss through this evaporation is dQE = -LMdQr Substitution of the individual heat sources and sinks into the first law of thermodynamics and dividing by Mdt yields dqs 1dM dQ dT -L- - [L(q - qe) + C (T - T )] - dM L = C dTdt e M dt pdt dt T97 (9) d v dp p dt 97

The saturation mixing ratio is given by q= - p ~ SE = 0.622 -s p" - e and the saturation vapor pressure by Teten's formula a(T 273.16) eS= 6.1078 exp a (T - 2b) (mb) (Murray, 1967) with a = 17.27 and b = 35.86. Differentiation with respect to time and combination of these two equations yields dqs qs (~ + qs) a(273.16 - b) dT dt (T - b)2 and TdT dv dp dp p dt dt substituting into Eq. (9) gives gw + [L(qs - qe) + Cp(T - Te)] + L dT dT - rq (c + q ) + a(273.16 - b) (10) L _, ~+ C (T - b)2 P THE ENVIRONMENT i) Environmctental air temperature Te: The environmental air temperature is assumed to decrease with constant values of lapse ratey e, under this assumption. Temay be expressed in the form T = T - y z. e eo -e 98

ii) Environmental vapor mixing ratio q: The saturation mixing ratio qs in the environment can s e be obtained from Teten's formula dq qse ( + qse) a(273.16 - b) sdze se se dz -'s(T - b)2 If the relative humidity RH is assumed constant, then q = RHqs, it follows that dq qe( + 1/RH q) a(273.16 -b) dz e 2 dz - le E: (T - b)2 Thus if the lapse rate yE, the temperature Teo and the relative humidity RH are initially chosen, T and q may H e e be obtained. RAINFALL RATE When the drops grow large enough to have an important terminal velocity (Vt > W), they fall relative to the updraft, the rainfall rate contributed by drops falling with terminal velocity at various altitudes is given by R = f D3N e V dD. 1 D6 o t,z c Integrating the equation, we obtain rrN P12 b_ m( R (- l) r m m!l (e XDbDm-r 1 6 P n r=lO (m-r)(-X)r+l b e D 99

where m = 3, 4, 5, 6 n = 0, 1, 2, 3; m = n + 3. 100

REFERENCES Austin, J.M., and A. Fleisher, 1948: A Thermodynamic analysis of Cumulus Convection. J. Meteor., 5, 240-243. Byers, H.R., 1965: Elements of Cloud Physics. The University of Chicago Press, 191 pp. Dingle, A.N.; 1970: Rain Scavenging Studies, Progress Report No. 7, Contract No. AT(11-1)-1407, U.S. Atomic Energy Commission, the University of Michigan, Office of Research Administration, Ann Arbor. Gunn, R., and G.D. Kinzer, 1949: The Terminal Velocity of Fall for Water Droplets in Stagnant Air. J. Meteor., 6, 243-248. Haltiner, G.T., 1959: On the Theory of Convective Currents. Tellus, 11, 4-15. Kessler, E., 1969: On the distribution and Continuity of Water Substance in Atmospheric Circulations. Meteor. Monogr., 10, No. 32, 84-88. Marshall, J.S., and W. McK. Palmer, 1948: The Distribution of Raindrops with Size. J. Meteor., 5, 165-166. Murray, F.W., 1967: On the computation of saturation vapor pressure. J. Appl. Meteor., 6, 203-204. Weinstein, A.I., 1970: A Numerical Model of Cumulus Dynamics and Micrographysics. J. Atmos. Scie., 27, 246-255.:101

PUBLICATIONS COO-1407-31. Trace Substances in Rain Water: Concentration Variations During Convective Rains and Their Interpretation, by D. F. Gatz and A. N. Dingle, Tellus, 23, 14-27, 1971. COO-1407-36. Comment on "Atmospheric Turbidity over the United States, 1961-66", J. Appl. Meteor., 10, 165, Feb., 1971. C00-1407-38. Rain Scavenging Studies, by A.N. Dingle. Progress Report no. 7, Contract no. AT(11-1)-1407, U.S. Atomic Energy Commission, The University of Michigan, Ann Arbor, COO-1407-38. 102

PERSONNEL 1. Mr. Yean Lee, graduate research assistant. 2. Mr. Duane Harding, graduate research assistant. 3. Mr. Bhanuprasad M. Joshi, analytical chemist, on deputation leave from Tata Institute of Fundamental Research Bombay, India. 4. Mr. James Fairobent, graduate assistant in research 103

UNIVERSITY OF MICHIGAN 3 9015 02841 2560