Mathematical Models Predict .Concentration-Time Profiles Resulting from Chemical Spill in a River W. Brock Neely" Health & Environmental Research, The Dow Chemical Co., Midland, Mich. 48640
Gary E. Blau, Turner Alfrey, Jr. Physical Research Laboratory, The Dow Chemical Co., Midland, Mich. 48640
With the increased use of the nation's waterways for the transportation of materials, there is an increase in the probability of spills. Once such a spill has occurred, there is an immediate need to predict the concentration profile of the chemical as the spill travels in order to assess the impact to both humans and the environment. This paper discusses the use of a mathematical model that has this predictive capability for common spills. Although the model is derived from the assumption that the chemicals are completely water soluble, it is also useful for partially soluble materials. The credibility of the model is demonstrated by comparing the concentration profile predicted with the actual profiles measured in two different incidents.
Chemical & Engineering Neus ( I ) published a report which predicted that the amount of chemicals shipped in this country by water will roughly triple by the year 2000. The frequency of accidents resulting in the discharge of the barge contents into water has been amazingly low. However, there is still the finite probability that these accidents will continue and result in the discharge of assorted materials into the nation's water. There are three common accidents that might occur which would result in the introduction of chemicals to a receiving body of water: (1) A barge could spring a major leak or buckle, thereby dumping the entire contents of the barge instantaneously into the river. This will be referred to as instantaneous loading of the chemical. (2) The leak could be small so that the chemical would enter the water a t a constant rate over a fixed interval of time. This is probably the most common accident and could refer either to a barge accident or to a leak from a point source located on shore. (3) There could be a combination of the above two examples. In such a situation an instantaneous loading might be followed by a slow infusion over a fixed time interval to the water or vice versa. Once such a spill has occurred or is occurring in the slow leak case there is an immediate need to predict the concentration profile of the chemical a t it travels down the river. The resulting concentration must then be matched with the known toxicological and other properties of the material so that appropriate action may be taken to alleviate any potential hazard. This report discusses a mathematical model that has this predictive capability for the three common types of spills. Although the model is derived from the assumption that the chemicals are completely water soluble, it is also useful for partially soluble materials. The credibility of the model is demonstrated by comparing the concentration profiles predicted with the actual profiles measured in two different spill incidents. Model Description T o a first approximation, a river may be visualized as a series of continuous stirred flow compartments as shown in 72
Environmental Science & Technology
Figure 1. In such a scheme, the output from each compartment is fed into the next compartment where the concentration of the output is the same as the concentration in the compartment. The rate of flow between compartments is related to the flow rate of the particular river in question. The minimum parameters needed to describe the river geometry by such a model are any combination of two of the following: cross-sectional area, flow rate in miledhour (velocity), and volumetric flow rate. A material balance for the flow of contaminant through the n t h compartment is given by the following differential equation:
where C, = uniform contaminant concentration in the nth compartment a t time t V , = volume of the nth compartment q = volumetric flow of the river assumed to be constant through each compartment he = rate constant for the evaporation of the contaminant in units of depthhime A = surface area of the compartment The evaporation of the material from the compartment is assumed to follow first-order kinetics, and the units are expressed in depthhime to make the equation dimensionally correct. With the availability of packaged numerical solution programs for solving differential equations of the type shown in Equation 1, it is not necessary to have a closed form of the equation. However, for the simple case of an instantaneous loading of a soluble chemical to a river, there are some useful relationships that can be derived from the solution of Equation 1, and equally important they can be handled with the aid of a simple calculator. In the first Compartment Cl(0) = M/V1 where M is the number of pounds of soluble chemical added instantaneously a t time t = 0. For all compartments other than the first, the concentration of contaminant initially is a t zero so that the initial conditions are simply: C,(O)=Oforn12
(2)
Dividing Equation 1 by V,, defining 8 as q/V,, and solving the resultant differential equation system with the
V = Volume Of Compartment w = Width I = Mixing Length q = Volumetric Flow Rate
Figure 1. Model of a river visualized as a series of continuous stirred compartments V , volume of compartment: w, width, 1, mixing lengths: q. volumetric flow rate
above initial conditions give the following expression for the concentration in the n t h compartment a t time t :
M (&),-I C , ( t ) = -exp -[(k,lh) V ( n - I)!
+ 8) t ]
(3)
This useful relationship demonstrates that by knowing the loading, flow rate, and evaporation constant, a concentration-time profile a t any position downstream may be readily calculated with a hand calculator. The one adjustable parameter is the mixing length. The mixing length is defined as that length of river in which the concentration of the chemical may be considered uniform. Future effort should be directed to develop a relationship between the mixing length and the river hydraulics. In the present analysis, the parameter was adjusted to obtain the best match with the experimental data. There are other relationships that may be derived from Equation 3: (1) The time for the maximum concentration to reach any point downstream is obtained by setting the derivative of the right-hand side of Equation 3 to zero and solving for t .
t,,, (2) Substituting t,,, mum concentration:
= ( n - l)/(ke/h + 8)
for t in Equation 3 gives the maxi1
8
(4)
27r ( n - 1)
Table I . River Geometry for Soldier Creek Used t o Make Concentration Predictions Miles from Velocity, Site source mPh Width, ft 2 2A 3 4
0.6 1 0.29 0.25 0.23
0.2 1.9 3.7 4.8
19.6 12.0 27.0 9.6
Notes: (a) data taken from Reference 2, ( b ) 2 Ib of soiuble dye added, (c) volumetric flow rate = 15.9 cfs, ( d ) dye was completely soluble a n d no evaporation was considered, a n d ( e ) mixing length estimated to be 60 ft.
32 Actual Data
Table I I . Geometry o f the Mississippi River Volumetric f l o w Velocity Width Depth
268,000 cfs 1.26 mph
4000 ft 36.3 f t
The other two types of accidents are much more difficult to handle in the closed form and do not yield the same simple relationships as discussed above. However, the resulting differential equations may be solved numerically on the computer. In the present case, as many as 121 simultaneous differential equations representing 121 compartments were successfully integrated in less than one minute on the IBM 370/158 computer using the IBM Continuous System Modeling Program (CSMP).
Testing the Model
Addition of a Soluble Chemical to a River. Bath e t al. (2) published an account of some studies they made on a small stream in Kansas. They added a soluble dye into the stream and measured the concentration-time profile downstream a t various points. This represents the situation for which Equation 3 is valid, and the river geometry extracted from their paper is shown in Table I. A mixing length of 60 f t fitted the data adequately. A comparison of the predicted concentration-time profile with the observed data is shown in Figure 2 where the experimental points were taken from Bath e t al. (2). The close agreement between the observed and predicted values is readily apparent, lending credence to the model for this situation. Addition of a Partially Soluble Chemical. On Sunday, August 19, 1973, a barge carrying three tanks of chloroform for midwestern terminals was damaged on the Mississippi River a t Baton Rouge, La. While it was being repaired, the contents of two 70,000-gal tanks were lost. T h e first tank ruptured a t 2:40 p.m. releasing its entire contents in a short period of time. The second tank sprang a leak a t 1O:OO p.m. on the same day and released its contents over a 45-min period. A total of 1.75M lb of chloroform was lost. The Louisiana Division of The Dow Chemical Co. USA began sampling the river a t 16.3 mi and 121 mi (New Orleans) from the point of spill. They determined the shape of the wave, the maximum concentration observed, and the time for the peak to arrive a t the two indicated points. In addition they determined that the chloroform was evenly distributed across the river a t a point 17 mi from the spill. The river geometry was taken a t the time of the accident and is shown in Table 11. These data compare favorably with information published by Everett ( 3 )on the hydrologic condition of the Mississippi River. In considering a volatile agent such as chloroform, an estimation must be made of the evaporation rate constant. For this purpose the data of Dilling e t al. ( 4 ) will be used. Briefly, these investigators measured the evaporation of several chlorinated solvents from a 250-ml beaker. The solution height was 2.55 in. The loss of chloroform from this experiment is shown in Figure 3. Experimental points taken from this figure were fitted by Equation 5 giving a rate constant of 0.364 ft/hr: In Ct = In Co - ( k / h ) t
(5)
HC”,5
Figure 2. Concentration-time profile of a soluble dye added to Soldier Creek in Kansas Continuous line represents prediction from model
Since this value represents evaporation from pure water and the Mississippi contains many things besides water, a modified value of k was desired. Again, Dilling et al. ( 4 ) atVolume 10, Number 1, January 1976 73
0.8
---
0.7
Pure Water
In The Presence Of 500 ppm Peat MOSS
0.6
0.5
G
0.4
0 0
0.3
0.2
0.1
0
0
20
40
60
80
100
Minuter
Figure 3. Evaporation of chloroform in water and contaminated water (Reference 4)
tempted to represent contaminated water by measuring the rate of evaporation in the presence of several possible contaminants. A value calculated in the above manner for the evaporation from a solution containing 500 pprn of peat moss gave a k in ft/hr of 0.255. The report of Everett (3)indicates that 300-400 ppm of sediment is a reasonable value for all flow conditions in the Mississippi. The author also noted that this would be higher for low flow conditions and since the chloroform spill did occur under low flow, a value for the evaporation of CHC13 of 0.255 ft/hr will be used in the modeling work. In addition to the evaporation rate of the agent, other properties that become important are such items as the water solubility, density, and partition coefficient. In the case of chloroform the water solubility is OB%, and the chloroform has a density greater than water. The partition coefficient between n-octanol and water is only 93 ( 5 ) .This value would agree with the observation of Dilling e t al. ( 4 ) that chloroform did not seem to be readily absorbed to organic matter. A t present, we are not able to make any ab initio calculations as to how these properties are related to what percent of the chloroform added to the river remains in solution and what percent remains as an insoluble mass on the bottom. Consequently, two approaches were made in predicting concentration profiles: (1) An assumption of complete water solubility will give predictions which represent a “worst possible situation.” The results of modeling this type of situation are shown in Figures 4 and 5. From these figures it is noted that the maximum predicted concentration a t 16.3 mi from the spill is 1.95 ppm, while a t New Orleans (121 miles from the spill), the maximum concentration was 0.625 ppm. ( 2 ) By introducing additional mechanisms an attempt was made to arrive at a model which would match the observed data as closely as possible and hopefully provide some plausible explanation for the observed phenomenon. In this approach a four-step iterative model-building procedure was used: (a) A sound physicochemical-hydrodynamic explanation was postulated to describe the manner by which chloroform entered the river flow system. (b) Differential material balance equations for the first compartment were written to describe the physical phenomenon. This involved subdividing the first compartment into sections which will be called layers for clarity. 74
Environmental Science 8, Technology
(c) The equations were incorporated into a CSMP computer program. Then, the constants or “parameters” of the model were adjusted in an attempt to make the model equations predict the concentration-time data collected a t the Plaquemine site 16.3 mi downstream from the spill and a t New Orleans 121 mi downstream. (d) If a lack of fit existed between the predicted and observed values, a new mechanism was postulated to accommodate these inconsistencies and the procedure returned to step (b). In following this model-building procedure, the principle of parsimony was adhered to rigidly. That is, the postulated mechanism and associated models were kept as simple as possible p n d then gradually made more complex until a further increase in complexity was not warranted by the data. Such a parsimonious procedure will result in a relatively simple model but does not preclude the existence of other suitable models. That is, the model generated is not necessarily unique. It is valid only in that it explains the concentration time data generated. As such, one must be very careful in trying to use the model for predictive purposes in any situation other than for the special Mississippi River spill situation for which it was developed. Although it would be an interesting exercise in model building, for the sake of brevity, only the final model developed will be presented here. For more details on model building see Reference 6. ---” wb 1800
-
1600
-
17
1200 1400
P
8
-
y
1000
-
800
-
SOP
-
Actual Data
U
200 400
0
. 10
20
30
40
5c
60
Elapsed Time, Hours
Figure 4.
Predicted concentration time profile assuming total solub,ili-
ty (16.3 mi from spill)
8oo 700
c 0
80
100
120
14C
160
180
Elapsed Time. Hours
Figure 5.
Predicted concentration time profile assuming total solubili-
ty (121 mi from spill)
350
3w Flow of river
I
1
I k,
-
Actual Data
- 250
Boundry layer
-s
Bottom layer of Chloroform
0
Figure 6. Hypothetical reactions taking place in the 1st compartment during the chloroform spill
200
0
2 150
u u
The spillage of chloroform in the Mississippi River is best understood by considering three distinct time phases: (1)The instantaneous ruptures of the first tank to the start of the infusion of chloroform resulting from the rupture of t h e second tank. This was a period lasting 7.33 hr. (2) During the infusion of chloroform from the second container; a period lasting 45 min. (3) After the infusion from the second container has been completed. Assume that all the chloroform holdup occurs in the vicinity of the spill-Le., in the first compartment of the discretized model system. This conclusion is supported by the analytical work performed by the Louisiana Division of The Dow Chemical Co. USA. Their analysis of the amount of chloroform that passed Plaquemine (16.3 mi downstream from the spill) and New Orleans, (121 mi from the spill) indicated that the difference could be accounted for by evaporation. The postulated model also precludes the existence of a dynamic equilibrium between the chloroform and the mud in the river bottom, a t least for this portion of the river. This is supported by the observations of Dilling e t al. ( 4 ) who claimed little or no binding of CHCls to sediment. Suppose that some fraction f l of the contents, L1, of the first tank remains as chloroform and, because of the high density of chloroform, drops to the bottom of the river as soon as the first rupture takes place. The other fraction 1 - f l dissolves uniformly throughout the rest of the first compartment giving some initial concentration in the river of L1 (1 - f l ) / V , , where V1 is the volume of the first compartment. I t is this amount which forms the wave front which proceeds downstream. The mass of chloroform from the spill ( L l f l ) on the bottom of the river, slowly diffuses with a rate constant k b into an adjoining boundary layer. From there the chloroform is transported into the river with a rate constant kp. These layers are depicted schematically in Figure 6. During the second stage of the accident the rate of infusion is:
ko = Ldtin where tin is the time period during which the infusion took place and L2 represents the contents of second tank. In this phase (1 - f z ) represents the fraction of ko which directly enters the water, and f z the fraction which enters the bottom layer. The differential equation related to each phase of the accident will now be described: (a) Phase 1. The material balance for each layer is shown below: For the bulk flow layer
For the boundary layer
100
50
I
I
I
I
where C1, Cp = chloroform concentrations in the bulk and boundary layers respectively V I , V2 = volumes of the bulk and boundary layer = rate constant for transport of the chloroform k2 from the boundary layer into the main compartment kb = rate constant for diffusion of chloroform from the bottom layer to the boundary layer B = Amount of chloroform in the bottom layer R1, k,, and h have been defined previously (b) Phase 2. During this time period, the material balance for the three layers becomes:
(c) Phase 3. The amount of chloroform in the bottom layer B2 a t the start of this period can be calculated from the integrated form of Equation 11: f2ko kb
+
BZ= - ( B 1 where Bi
-e)
exp [-kb(t
- 7.3)]
(12)
= amount of chloroform in bottom layer a t the
start of phase 2 = time a t which phase 2 begins The differential equations describing the chemical flow during Phase 3 are the same as Equations 6 and 8 with the single exception that Bo is replaced by Bz. The equations developed above for modeling the first compartment were simulated on the IBM 370/155 computer using CSMP. The parameters were adjusted until the concentration-time profiles shown in Figures 7 and 8 were generated corresponding t o the data collected a t mile 16.3 and mile 121 respectively. The close agreement between the calculated and observed values lends credence to the model representation. The following parameters produced the two plots: (t
- 7.3)
(7) f l = 0.82 kp = 1.0 hr-'
For the bottom layer
dB -= dt
-kbB
f p
= 0.97
kb
= 0.003 hr-l
The model is extremely sensitive to the values of f l and f 2 . However, there is significant interaction between the rate Volume 10, Number 1, January 1976 75
c
300 I
P
250
:200
.
I
*Actual Data
.
.*
The above model fits the data adequately and provides a plausible explanation of the chloroform spill in the Mississippi River. It demonstrates how the soluble contaminant concepts can be modified via an interactive model-building procedure to arrive a t a model suitable for describing a contaminant with different physicochemical properties. Future work should be directed a t estimating these parameters from the physical chemical properties of the spilled chemical.
Acknowledgment 50
60
80
~
100
’
120 Elapsed Time, Hours
140
160
’
180
Figure 8. Concentration time profile of chloroform (121 mi from point of addition)
constants k2 and kb such that their ratio is quite sensitive although their absolute values may change as much as 25%. The values of f l and f 2 indicate that the majority of the chloroform goes into the bottom layer. In particular, 97% of the contents of the second tank pass into this layer even though the rupture took place over a longer time period. The absolute values of k z and kb are not too important since they are so intimately involved with the particular environment of the accident. However, the ratio k2/kb = 3000 indicates that the diffusion of chloroform into the boundary layer is considerably slower than the transport from the boundary layer into the bulk flow layer.
The analysis by G. W. Daigre and E. J. Brown on the chloroform spill was most complete. These people and the many others in the Louisiana Division of The Dow Chemical Co., USA, are to be congratulated for an excellent task. Without these data and the many helpful discussions with Daigre, this study would never have been possible.
~
~
~
References (1) Chem. Eng. News, p 7, )feb. 25, 1974). (2) Bath, T. D., Vandegrift, A. E., Hermann, T . S. J . Water Pollut. Control Fed., 42,582 (1970). (3) Everett, D. E., “Hydrologic and Quality Characteristics of the Lower Mississippi River,” published by the Louisiana Department of Public Works, Baton Rouge, La., 1971. (4) Dilling, W. L., Tefetiller, N. B., Kallos, G. J., Enuiron. Sci. Technol., 9,833 (1975). (5) Leon, A., Hansch, C., Elkins, D., Chem. Reu., 71,525 (1971). (6) Blau, G. E., Neely, U’.B., Adu. Ecol. Res., 9, 133 (1975).
Received for review March 28, 1975. Accepted September 18, 1975.
Monitoring California’s Aerosols by Size and Elemental Composition Robert G. Flocchini,” Thomas A. Cahill, Danny J. Shadoan, Sandra J. Lange, Robert A. Eldred, Patrick J. Feeney, and Gordon W. Wolfe Crocker Nuclear Laboratory and the Department of Physics, University of California, Davis, Calif. 95616
Dean C. Simmeroth and Jack K. Suder California Air Resources Board, 1709 1 1th St., Sacramento, Calif. 958 14
The atmospheric aerosol consists of a complex ensemble of particles in an infinite combination of physical and chemical states. Despite their importance in reducing visibility. affecting human health, and soiling materials, their complexity has hindered attempts to include detailed information on aerosols in air-quality monitoring programs. Generally, only the total suspended particulate present a t a site during a 24-hr period is measured. Some information on chemical composition is extracted from aerosol samples, but analytical costs limit such analyses to a few important species on representative samples. Recent advances in energy-dispersive X-ray analysis have resulted in dramatically reduced costs for quantitative, multielement analyses of air samples. One can therefore visualize more complete aerosol monitoring efforts, including information on particle size and elemental content of aerosols at many locations for extended time periods. Such a program has been established by the California Air Resources Board, working in conjunction with the University of California, Davis. Up to 15 sites were selected a t locations representative of large areas of the state. Aerosol samples were collected in three particle size ranges by 76
EnvironmentalScience & Technology
means of Environmental Research Corp. Multiday Impactors. These units are rotating drum impactors of the Lundgren type with after filters ( I ) . Once a week samples were sent to Davis and analyzed by ion-excited X-ray emission for sodium and heavier elements. This paper will describe the analytical methods used in the collection and analysis of the aerosol samples, with emphasis on validations used to ensure accuracy.
Site Selection The sites were selected in such a way as to present aerosols typical of larger areas throughout California. The sites are shown in Figure 1 and, from north to south, they are: Geyserville, Sacramento, Richmond, Oakland, Livermore, San Jose, Salinas, Bakersfield, Azusa, Los Angeles, Riverside, Los Alamitos, Indio, and El Cajon. An ARB Mobile Air Surveillance Unit was also equipped with one of the samplers. Sampling Methods Aerosol samples were collected by means of Environmental Research Corp. Multiday Impactors (presently manu-
~