Environ. Sci. Technol. 1988, 22, 1325-1336
Complexing Equilibria in Mixed Ligand Systems: Tests of Theory with Computer Simulationst Donald S. Gamble” Land Resource Research Centre, Agriculture Canada, Ottawa, Ontario K 1A OC6, and Department of Chemistry, Concordia University, 1455 Boulevard de Maisonneuve West, Montreal, Quebec H3G 1M8
Cooper H. Langford Department of Chemistry, Concordia University, 1455 Boulevard de Maisonneuve West, Montreal, Quebec H3G 1M8
The equilibrium complexing of an ion or molecule by a mixture of binding sites is a general type of chemical system that has widespread geochemical and biological importance. The analytical chemistry of such systems must include methods for determining the chemical forms of the small ions and molecules, in addition to their total amounts or concentrations. Two parallel streams of closely related theoretical work have developed independently. One body of theory has developed in the geochemistry of soils and natural waters for cation reactions with humic materials. The other is found in biochemistry, with peptides commonly providing the mixture of binding sites. A number of misconceptions about the theories still persist. In this paper, computer simulations are reported for ligand mixtures for testing and comparing the theories. Some experimental methods are proposed for avoiding the practical limitations of the calculation techniques.
Introduction The equilibrium complexing of ions or small molecules by a mixture of ligand sites confined to locations on polyelectrolyte molecules represents a generic type of system that is important because of its biological ubiquity. The central characteristic of this generic type of equilibrium system is that a given reagent is complexed by a mixture of chemically nonidentical ligand sites. The significance of the generic type of complexing system is indicated by the fact that it has been independently investigated in two unrelated fields of research. One of these is the geochemistry of soils and natural waters, in which humic polyelectrolytes provide the complexing sites. The second is biochemistry, with the polyelectrolyte quite commonly being a peptide. A number of theories and models have been developed for calculating and interpreting these mixed system equilibria, and various reviews are available (1-5). Sposito has published the most recent review (6). It is an excellent one because it uses underlying concepts with which to give a unified discussion of the various theoretical mathematical models and theoretical mathematical descriptions. In all of the cases covered, the law of mass action has been applied macroscopically to whole mixtures. Marinsky, Gupta, and Schindler also followed this historically conventional practice in their investigations of the polyelectrolyte chemistry of humic-metal ion reaction equilibria (7,8). In 1926, for example, Simms proposed that the law of mass action be applied to the overall equilibria of a whole mixture (9). He called the resulting equilibrium function a “titration constant”. From the earliest years of metal ion-fulvic acid research, it has been common practice for investigators to calculate and present their experimental
* Author to whom correspondence should be addressed at the Land Resource Research Centre. Contribution No. 87-69. 0013-936X/88/0922-1325$0 1.50/0
results by applying the law of mass action in just this way to the whole fulvic acid mixture. The early authors were, however, unaware of Simms’ original proposal (10-20). Three important points should be noted about this historically customary practice. The first is that during the early years of the research, no theoretical or practical alternatives were available. Second, it was universal textbook convention to think about the law of mass action in terms of mathematical constants, and early workers were made uncomfortable by the experimentally obvious fact that the complexing of metal ions by fulvic acid did not produce constants. The final point is that this historical practice does, as Simms recognized, have a great deal of practical merit. A clear-cut requirement follows logically from this. The nature of the equilibrium functions produced by the Simms procedure must be properly understood. The key insight is the concept that the application of the law of mass action to mixed equilibrium systems produces mathematical functions, rather than mathematical constants. The theories and models should be considered from this point of view. Klotz and Hunston referred to two types of models (21). One consists of empirical mathematical functions used strictly for statistical correlations of experimental data. The other type uses a mathematical formulation based on molecular level postulates, which have been adopted for subsequent experimental testing. Perdue (5) used this type of model for developing better insight a t the molecular level. There has been an important type of activity that is complementary to modeling. The present work deals with the two streams of research that have developed in geochemistry and biochemistry on rigorous theoretical analysis. Although they developed independently, they have some important features in common. (a) Experimental parameters such as site coverage and equilibrium functions measured on whole mixed samples have been expressed as weighted averages. (b) Various mathematical methods have been devised for extracting equilibrium constants from the experimental weighted average functions: differential equilibrium function calculations; Shuman’s Laplace transform method (22);Klotz and Scatchard graphical plot methods (23,24); Sposito’s quasiparticle model calculations (6). Klotz and Hunston discussed limitations of the graphical plot methods (24). They indicated, for example, that the use of Scatchard plots for ligand mixtures may yield parameters that are complicated combinations of component equilibrium constants instead of equilibrium functions for the whole system or for individual components. Sposito has reexamined this theoretically however, in terms of the Scatchard quasiparticle model analysis (6). Differential equilibrium functions, Laplace transform methods, and quasiparticle model calculations deserve further examination.
0 1988 American Chemical Society
Envlron. Sci. Technol., Vol. 22, No. 11, 1988
1325
The formal recognition that metal ion-fulvic acid complexing equilibria are characterized by mathematical functions rather than by mathematical constants was the factor that led to the weighted average equilibrium function theory (25). Its original purpose was to explore and elucidate the inherent chemical meaning of the type of equilibrium function that was already in general use throughout the literature. The important point is that the weighted average equilibrium function theory was not the invention of a new type of treatment for equilibrium data. It simply sought to explain one already in general use. In that sense it is conceptually conservative. The basic features of the theories and calculation methods have not always been well expressed, and confusion has arisen on some points. They include the following topics: (1)It has been widely assumed that simple graphical plot methods similar to Scatchard plots may be used for recovering the equilibrium constants of components from the experimental measurements on mixtures. The linear plot methods used by Klotz and Scatchard were intended for the calculation of a single equilibrium constant, representing a monomeric pure reagent, from slopes and intercepts. Mixtures, especially of polyelectrolytes, tend to give nonlinear plots. This has led some authors to make arbitrary decisions about the variable slopes and about the types of numbers of empirical parameters with which to represent the system. Perdue has pointed out that the parameters obtained in this way cannot be assigned any real chemical significance ( I ) . Klotz and Hunston found (24, 26) that some of the parameters obtained are equal to complicated combinations of individual equilibrium constants. Sposito (6) emphasized that practical methods for calculation and interpretation must account for these facts. The number of unknowns is usually greater than the number of independent equations. For this reason, tests for the recovery of Ki values have been confined to the weighted average and Laplace transform methods. (2) The weighted average equilibrium function, R, defined and calculated for mixed equilibrium systems, has frequently been regarded as “a constant that failed”, It has been shown instead that it is inherently a mathematical function rather than a mathematical constant. The independent variables include especially the chemical composition variables of the mixed system. It should be noted that none of the other theoretical models or descriptions account rigorously for all of the chemical composition variables. Some recent computer models have been based on just one or two key empirical or semiempirical parameters. While experimental testing of these may reveal some useful approximations for particular cases, the need for generality still exists. The complete list of chemical composition variables employed in a rigorous total system description may be the best means of achieving this. Finally, none of the other theories can be related as rigorously to conventional Lewis and Randall chemical thermodynamics through these variables. (3) It has recently been suggested that the differential equilibrium function does not properly represent individual chemical species, and therefore its use is incorrect. Also, some authors still consider it improper for a law of mass action equilibrium function to vary with experimental conditions. There are three responses to these points of view. The first concerns the original purpose of the weighted average equilibrium function theory. In geochemical systems outside living organisms, the biologically important variable is often only the thermodynamic activity of the small ion or molecule in free form. In that 1326
Environ. Sci. Technol., Vol. 22, No. 11, 1988
case, the geochemical ligand mixture is only relevant as a complexing system controlling that activity. Individual species in the ligand mixture are not separately important to the biological problem. The weighted average equilibrium function theory was therefore developed as a description of mixed ligand systems, for the purpose of predicting the concentrations and thermodynamic activities of free small cations. There was neither practical need nor any intention to use it for the identification of complex species. For this reason, the original presentation of the weighted average equilibrium function theory did not refer to species. Instead it was formulated (25)in terms of “.., the ith infinitesimal increment of carboxyl group molalities at a given value of c; = mAt + mA;H which is equilibrated between ionized and unionized forms ...”. The question about whether or not the differential equilibrium function does identify species in the mixture therefore relates to an objective outside the original scope of the work. Complex species might, however, be important in biochemical systems. The response for objections to variable equilibrium functions is that their routine use has been accepted in the literature for many years. The earliest investigators published equilibrium functions for humic-metal ion complexing. In addition, AGO and AHo have always been presented as functions of sorption site coverage in the surface chemistry literature. In fact, both differential and integral functions are used for site coverage. The third point is that the original scope of the research can now be extended to an investigation of the relationship of the differential equilibrium function to the Kiequilibrium constants of components. This can be done with computer simulations. Perdue has demonstrated that computer simulations may be used for investigating the Cu(I1)-complexing properties of mixed ligand samples ( I ) . Although his objectives were different from those above, his strategy can be profitably adapted to the present purposes. (4) It should not be assumed that the variation of K with increased cation loading was caused by the electrostatic fields of polyelectrolyte molecules. The original weighted average equilibrium function theory postulated instead that such variation is an inherent property of mixed equilibrium systems (25). Any polyelectrolyte characteristics are superimposed as secondary effects. This can be tested with computer simulations. (5) The weighted average and differential equilibrium functions have sometimes been used as if they were interchangeable, especially for estimates of AGO. Their published experimental values are quite commonly different, and this property of mixtures can also be theoretically deduced and demonstrated with computer simulations. Rigorous, quantitative theory should be available both for guiding decisions about the types, amounts, and arrangements of analytical chemical measurements to be made and for use in the interpretation of the resulting analytical chemical data. Provided that activity coefficients are properly accounted for, then the calculation of AGO from K is strictly correct, and in accord with Lewis and Randall thermodynamics. The AGO therefore has the usual chemical meaning, but relates to a very small portion of a mixture. A AGO that has been calculated from K is ambiguous and does not represent any individual component complex. It is not a conventional function and is not in accord with Lewis and Randall thermodynamics. It could not be interpreted at the molecular level. For this reason the calculations in this paper are intended to serve three purposes. The first purpose was to identify the common characteristics of the two bodies of theory. The
second purpose was to correct the misconceptions. The final objective was to determine the practical limitations of calculation techniques.
4
i-1
cu2+
+ L 5 CUL
IT: =Mc/McUM~
(2)
The subscripts c, indicate components of the mixture. In the same way, eq 3 describes the total concentration of free ligand.
ML = ML1 + ML, + ML, + ML4
(3)
The material balances for total Cu(I1) and total ligands are accounted for by eq 4 and 5. The mole fraction of total ligands occupied by Cu(I1) is defined by eq 6 The mole
ccu = MCU + Mc CL = ML + Mc xc
= Mc/CL
= ML/CL
b
b
dissociation equilibrium of fulvic acid (25). mAiH is the molal concentration of the unionized portion of the ith component of the mixture. Subsequently, Klotz and Hunston independently proposed eq 16 for association
equilibria (26). y is an integer. Equation 16 represents a general category of averages that is a conventional topic in statistical theory. Within the general category, the integer defines particular types of averages. Typical values for y include 0, 1, 2, etc. Klotz and Hunston did not determine what value of y is valid for mixed ligand equilibria. For the particular case of y = 1,the equation of Klotz and Hunston becomes essentially the same as eq 14. niis (moles of ith type site/mole of protein molecules) that still remain free. ( h ) , is the particular type of average defined by a given value of y. Equation 17 describes the amount of free ith ligand as a mole fraction of the total amount of all the ligands.
(4) (5) (6)
fraction of total free ligands is defined in the same way by eq 7. Each ligand in the mixture comes to its own XL
(14)
(1)
The total molar concentration of complexed Cu(I1) is given by eq 2.
Mc = Mc, + Mc, + Mc, + Mc,
i=l
type was eq 15, which is the exact analogue for the H+
Theory I. Weighted Average Equilibrium Function for a Mixture of Four Ligands. The first case to be considered is that of equilibrium Cu(I1) complexing by a sample containing a mixture of four simple monomeric ligands. Following the historic convention for fulvic acid, eq 1 is obtained by applying the law of mass action to the whole sample all a t once.
4
K = CK;ML,/ C MLi
(AX& = ML~/CL
(17)
(Ax& is therefore the ith small portion of the total mixture of free ligands. It is explained by the material balance eq 18. The notation chosen here reflects the fact that in the n
(7)
equilibrium in its reaction with Cu(I1). The ith such case is described by eq 8. This ith case has its own material
limit
K
cu2+ + L, & CUL,
Ki= MC,/McUM~, (8)
balance equation for total ligand, as shown in eq 9. The
CLi = ML, + MC&
(9)
K defined and calculated for the whole system in eq 1can be related to the Ki values of the component ligands by using material balance relationships. Substituting eq 2 and 3 into eq 1 gives eq 10.
Equation 13 may now be converted to mole fraction form by using the relationships of eq 7 and 17. The result is the pair of eq 20. According to eq 20, K is not a constant,
+ K ~ ( A x L+) ~K ~ ( A x L+) ~K ~ ( A x L ) ~ K = K~(AxL)I + ( A x L ) ~+ (AX,), ( A x L ) ~+ (204
i K ~ ( A x L+) ~& ( A x L ) ~ + K ~ ( A x L ) ~ K = K i ( A x ~ )+ XL
(20b)
From eq 8 one gets the rearranged form, eq 11.
Mci = Ki Mcu MLi
(11)
Equation 12 is.then obtained by introducing the relationships of eq 11 into eq 10. Simplifying this gives eq 13.
K= KIMC~ML, + K ~ M C ~ M+L&Mc~ML, , + KMC~ML, MCu(ML, + ML, + ML, + MLJ (12)
K ~ M L+, KzML,+ K ~ M L+, K ~ M L ~ K= (13) ML1 + ML, + ML, + ML4 I t is often useful to write eq 13 in the condensed form of eq 14. The first law of mass action equation of this general
but instead is a weighted average of the K;values of the four component ligands. A plot of K vs XL produces an envelope function
K
= f(XL)
(21)
xL,the independent scan variable, is the mole fraction of the whole ligand mixture that still remains free, as defined by eq 7. It must now be determined whether or not the original K , values can be recovered from the envelope function, eq 21. This is done by examining in detail eq 22, which is obtained by rearranging eq 20.
K X L= K ~ ( A x L+) &~ ( A x L ) ~ + & ( A x L ) ~ + & ( A x L ) ~ (22) The compleximetric titration of the mixture with Cu(I1) decreases each free ligand in its turn to zero, with the strongest complex always forming first. Equations 23 and Environ. Sci. Technol., Vol. 22, No. 11, 1988
1327
( K d i = K ~ ( A x L+) ~K ~ ( A x L+) ~& ( A x L ) ~ (24) 24 thus apply to the strongest ligand. Equations 25-27
(Rxd = ( R x ~ ) ~
-.
lim
(AXLA
0
(25)
then become valid as the second ligand is all occupied by Cu(I1). Since (AXLA 5
(Ax,),
of component K, values from titration data is imposed by the resolution of a practical titration. Just as with K , a plot of KxL against xL gives a continuous curve. It is in fact another envelope function. This is expressed by eq 36. It is simply a rearranged form of the envelope function
(26)
(27) (RXL):! = K ~ ( A x L+) ~& ( A x L ) ~ An experimental plot of KxLvs xL produces another envelope function scanning the titration. According to eq 3, 7, and 17, XL may be related to component ligands by eq 28. If the first ligand is all reacted this reduces to eq = ( A x L ) ~+ (AX,), + (AxL)~+ (AX,), (28) 29. Likewise the removal of the first two ligands gives eq
-.
lim
AXL
0
(Ax,), + (AXLI3 + (AXL14
(29) 30. Equation 31 then indicates the region of the titration (XLh
=
@XL)3
+ (AXL)4
(30)
AXL = (XL), - (XLA
(31) scan between the removal of the first ligand and the removal of the first two ligands. The increment of the envelope function (RxL)over the scan region from ( A x L ) ~to may now be derived. Subtracting eq 24 and 27, = (AX,),
(Ex,), - (RXL), A(Rx)i,2 = K ~ ( A x L+) ~& , ( A x L ) ~+ & ( A x L ) ~ A(Rx,),,2 =
(32)
K3(AXL)3 - K3(AXL)4 (33) = K2(AXL)Z
(34)
A(RXL)l,2/(AXL)2
(35)
A(RXL)1,2
K2 =
A t this stage in the logical analysis, practical analytical chemistry intrudes on the pure mathematics. Usually two or more of the component ligands will react at the same time as the titration proceeds. The reason is that each component is separately distributed between free and complexed forms according to eq 9. Furthermore, all such component equilibria are linked to give a buffering effect. Consequently, (Ax,), might become quite small by the time (Ax,), vanishes. The farther the titration proceeds, the more difficult it becomes for ( A x L ) ~to remain smaller than (AxL),. This explains why eq 23 exists instead of a sharp cutoff. Also as the titration proceeds, the amount by which ( A x L ) ~remains smaller than (AxL), will become harder to detect above the resolution limits of the titration method. The closer the two K, values are, the more closely their (Ax,), values will tack together as the titration proceeds. This can increase the experimental difficulty of achieving the condition in eq 24 and 27. In an extreme case, if two components have their Kivalues too close together and are titrated too early in the titration, the points at which the conditions 24 and 27 are met may be off the front end of the titration curve. In effect, the first drop of titrant grossly overwhelms the very small difference in behavior of the two components. Further into the titration, the buffering effect of the ligands already titrated and those beyond will control the activity of free titrant sufficiently so that the small difference is no longer overwhelmed. A corresponding phenomenon will exist at the tail end of the titration, where there is insufficient buffering effect beyond. Another limitation in the recovery 1328
Environ. Sci. Technol., Vol. 22, No. 11, 1988
[A(~~xL)/AxLI,,=,;
= [ ~ ( K X L ) / ~ X L ~ ,(37) ,=~,~
K2 = [d(RXL)/dXLIXL=XL~
XL
(XLh =
(36)
K X L = f(XL)
eq 21. Provided that the conditions of eq 23 and 25 are met, eq 35 must define a slope that is a tangent at some point, XL = xL’on the envelope function between (XxL), and ( K x L ) ~K. z may then be recovered as the derivative by collapsing the small interval AxL onto this point, to give eq 37 and 38. The subscripts indicate that the slope is being calculated at the particular point xL = x,’.
(38)
11. Calculation of the Envelope Function from the Component K, Values. The computer simulation of a mixed ligand system consists of the calculation of thdt type of equilibrium function which would be determined experimentally by a compleximetric titration. This means that the weighted average equilibrium function K is calculated from the equilibrium constants of the components over a range of site coverages. According to theory, K may be quantitatively predicted from the chemical composition of the mixture and the component X, values. Cabaniss et al. (27) studied a number of theoretical models and descriptions from the points of view of chemical credibility and usefulness for practical calculations. The ultimate purpose is to interpret experimental R data. The inverse operation pursued here is the deduction of system K values from the K, values of known components. The calculation requires that the molarity of free Cu(II), Mcu, be expressed in terms of Ccu, the molarity of total Cu(I1) and the component K, values.
Ccu = Mcu (1 + K ~ M L+, &MI,,
+K ~ M+ L~ K~MLJ
(39) The necessary relationship is obtained by combining eq 4, 2, and 8. This leads to eq 40, which expresses the Mcu = Ccu/(l + K ~ M L+, K&hz
+ K ~ M L+, K ~ M L J
(40) concentration of free Cu2+in terms of the known composition of the sample. Equations 8 and 9 for the ith ligand are combined to give eq 41 and 42, ligand still in the free Mc, = CL, - ML,
K, =
(CL, -
ML~)/Mc~ML,
= C L , / ( ~ + K,Mcu) form.
(41) (42) (43)
For convenience, eq 44 is derived from this. (44) B, E K,ML,= K , c L , / ( ~+ K,McJ McU= Ccu/(l
+ B1 + Bz + B3 + BJ
Mc = ccu - Mcu ML = CL - Mc
(45) (46) (47)
Equation 45 is finally obtained from eq 40 and 44 for use in practical calculations. The concentrations of total complex and total free ligand are calculated from eq 46 and 44, which are simply the rearranged forms of eq 4 and 5. Mcu and the B, terms must be evaluated by iterative calculations before M, and ML can be calculated. This in turn permits the calculations of R, xc, and XL. Equation
44 is related to one of the calculation techniques described by Hunston (28). This method assumed a continuous spectrum n = n(K) of binding constants and employed a binding function. The binding function (moles of small molecule/mole of macromolecule) is a sum of terms each of which is closely related to the expression eq 44. Hunston’s calculation was, however, the inverse of the one described here. 111. Extension of the Theory to a Hypothetical Sample with 100 Ligands. For a hypothetical sample containing 100 ligands, eq 2,3,20,35, and 45 are extended as follows: 100
M, = EM,,
(48)
1=l
100
ML = EM,,
(49)
1=1
100
100
K = CKl(A~L)l/E(A~L)l 1=1
i=l
(50)
100
R = CK,(AXdi/XL
(51)
1=1
Kl =
A(RxL),-I,,/(AxL)~
(52)
100
Mcu = CCU/[l
+ CPlI i=l
(53)
A natural sample, for example of humic material, may have a very large number of similar components. Each is likely to be a very small fraction of the total sample. Equations 48-53 may be easily generalized and extended for such natural samples. The number 100 is then replaced by the large number n, and each weighting factor is very much smaller than 1. IV. Recovery of the Kl Values. A . Differential Equilibrium Function Method. The plots of K vs xc in Figures 1 and 2 illustrate the fact that the calculation of K from four discrete component Kl values produces a smooth continuous curve. That is, it is an envelope function. The ExLvs xL plot is consequently also an envelope function. For this reason, the first derivative gives the envelope function K = f(xL)from which the recovery of component K, values put into the calculation are shown for comparison. When they exist, the xLaddresses of the component Ki values are also indicated. It is reasonable to anticipate that the effectiveness with which the component Kl values can be recovered and the sizes of the calculation errors will depend on the percent contributions of the component ligands to K. The percent contribution of each ligand to K is calculated according to eq 54, which is derived from eq 20. Y,= IOO-K,(Ax,), / K x L % (54) Envelope functions of Y,vs x, have been plotted for four ligands in Figure 3. The relationship actually observed empirically is that the two ligands for which Kl values are recovered have maxima in their percent contribution curves. The two not recovered do not. B. Laplace Transform Method. An alternative starting point exploited by Klotz and Hunston (24) for a description of metal titration of proteins and by Shuman and his students (22) for titration of humics begins from the rearrangement of eq 8 that uses the definition of xc = Mc/cL. This is xC= m c u / ( l + RMcu) (55) If K were a constant in eq 55 rather than a function, this expression could be used in the Shuman procedure. Since
400K1
-
30.0
7
P
,$ 20 0 -
we deal with a mixture, we can write for m components labeled by i. m XC
= CAxiKiMcu/(l i=l
+ KiMcu)
(56)
because (57)
Now Hunston and Shuman take a daring leap by considering a continuous distribution of ligands, a step that will be “effectively reversed” when numerical approximations Environ. Sci. Technol., Vol. 22, No. 11, 1988
1329
100.0
1
90.0
i=l
“““14 40.0
-
30.0
-
20.0
10.0-
,
I
0.10
I
0.20 0.30 0.40 0.50
I
I
0.60
0.70
0.80 0.90
I
1.00
xc I
I
I
0.90 0.80 0.70
I
I
0.60 0.50
I
0.40
I
I
I
0.30 0.20 0.10
1
0.00
ion complexing stability constants. A mathematical formalism may be adopted in the form of a collection of hypothetical complexing sites. This collection has the same number of sites as the real one. The real and hypothetical collections are both assigned the same hypothetical stoichiometry of a 1:l mole ratio of ligand sites to metal ions. The hypothetical collection has a single complexing stability constant that is equal to the average stability constant of the real one. Simple mathematical analysis has shown this to be a weighted average. A quasiparticle is a hypothetical polyelectrolyte molecule upon which there is one of the hypothetical sites. Sposito has emphasized three points. The first is that the Scatchard quasiparticle model is used because it may be conveniently applied to practical calculations for particular experiments. The second point is that it produces average stability constants. Finally, the quasiparticles and hypothetical complexing sites should not be mistaken for real chemical entities, and the resulting average equilibrium function calculated for the whole collection of nonidentical sites cannot be identified with any one real site. V. Weighted Average xo: Mole Fraction of Total Sites Occupied.
XL
Figure 3. -Mole percent contributions of components to the weighted average K . The effect of site coverage is shown.
Let 8, = Mct/CL,
(61)
are introduced to deal with real data. On this basis, they write
This is the mole fraction of the ith component that is occupied by Cu(I1). Correspondingly, eq 62 gives the mole (1 - 8,) = ML,/CL, (62)
and
fraction of this component that is free. The experimentally observed mole fraction xc may be expressed in terms of the mole fractions of the components occupied by Cu(I1). The substitution of eq 6 gives eq 63. The mole fractions xc
The integral equation cannot be solved in closed form, but a numerical approximation developed by Ferry (29) can be used to estimate Ax(K) given an experimental xc vs McU function.
= (MCI+ M,, + Mc,
+ Mc,)/CL
(63)
8, are introduced by using eq 61. The relationship to the components is then explained by the resulting eq 64. For x c = (llCL1 + 82CL2 + 83CL3 + 84CLJ/CL 1
4
= -CB,CL~ (64) cLL=1
Ax(K) =
the ith component, K, can be described in terms of the mole fraction 8,, by combining eq 61, 64, and 8 to obtain eq 65. Equation 65 is rearranged to give the mole fraction
K , = o,/Mcu(l - 8,) (65) in terms of the equilibrium constant, in eq 66. This gives 6, = McuK,/(l + McuK,) (66) where log a is chosen as 0.2. This system is solved with a computer program that first smoothes the experimental binding curves and then uses xc values at McUin a procedure to give Ax(K) vs log K , called an “affinity spectrum”. The ability of the numerical approximation to handle discrete components rather than a continuum was demonstrated by Shuman et al. (22) who used the equations to describe the components of the titration of H,PO,, HPO;-, and with protons. The area under each peak in the spectrum represents a distribution of components close in Ki.For a future set, the peaks are simply artifically broadened descriptions of the components but the quantity of each component is proportional to the area under the peak. C. Sposito’s Scatchard Quasiparticle Model (6). In this model, the whole mixture of complexing sites is conceptually divided into discrete categories. Each category is actually a collection of chemically similar, nonidentical complexing sites. They therefore have a range of metal 1330
Environ. Sci. Technol., Vol. 22, No. 11, 1988
eq 67 as an alternate form of eq 64. A simple extension 1 4
xC= -c1=1 C [McUK,/ (1 + McuK,)1CL,
(67)
to 100 component ligands gives eq 68. In 1971, Klotz and
Hunston presented a form of eq 67, shown here as eq 69. r noXc = ?[kiA/(l i=l
+ kiA)]ni
(69)
n, and no are (moles/mole of protein) of ith and total complexing sites. A is the molarity of the free small molecules or ions. ki is the complexing stability constant
Table I. Cu(I1) Complexing Stability Constants Used for Computer Simulations: Mixtures of Four Components Cu2+ + L & CuL
K = M,/(McuM~) 4
M
total concentration of all ligands: CL = CCLi= 5.0 X i=l
mixture component 1
CL (free + complexed), M
Li
1
malonate
2
carbohydrazide glycine ethyl ester tartrate
3
4
1.60 X 1.10x 10-3
0.80 X 1.50 X
Ki
3.548 X 8.318 x 1.380 X 1.585 X
lo6 104
lo4 lo3
for the ith type of sites. Klotz and Hunston were able to define the parameter r because they were dealing with identical protein molecules. In contrast to this, it is necessary to use xc for fulvic acid because no two polymer molecules in a sample are identical. In 1972, Gamble (25) independently deduced the exactly analogous eq 70 for the H+dissociation equilibrium of fulvic acid. 1 QA = --cpici CA 1
k
(70) + KAi)lci CAi=O Calculation of the I?. Envelope Functions. The purpose is to calculate the envelope functions K = f ( x J and K = f(xL)from the Kivalues of the component ligands and the chemical compositions of t h esample solutions. The independent chemical composition variables xc and xL must also be calculated. Table I lists the four simple monomeric Cu(1I) complexing ligands chosen from the literature. They have been chosen so that the four stability constants Ki span a range of -2 orders of magnitude. Ki is defined as shown at the top of the table, without proton competition. The choice of ligands was made from the tables compiled by Ringbom (31) and those of Martell and Smith (32). The quantities known from the sample compositions are as follows: (1)the stability constants Kiof the component ligands; (2) CL,the molarity of the total ligands; (3) CL,, the molarity of the total ith component ligand; (4) Ccu,the molarity of total Cu(I1). Step 1: Iterative Calculation of the Pi Terms and Mcu. Equations 44 and 45 were used for the iterative calculations, with zero-order approximation being obtained from eq 71. The iterations continue until successive McU = -x[KAi/(mH
Pi
Kic~,/(1+ KgCcu)
(71)
values differ by less than 0.1%. This might require 10 to 100 iterations. Step 2: Calculation of the Concentrations of Total Complex and Total Free Ligand. Once Mcu values were available, the concentrations of total complex and total free ligands, M, and ML, were calculated from eq 46 and 47. K,x,,and xL were then calculated from eq 1, 6, and 7,
Results a n d Discussion At least two properties of K can be anticipated in advance. Since a compleximetric titration must always form the strongest complexes first, it can be predicted that K will decrease as the mole fraction of complexing sites-occupied increases. It is also expected that a plot of K vs xc will produce a continuous curve or envelope function, even though K is defined in terms of discrete Kivalues. The curve in Figure 1confirms both of these expectations. It was calculated according to eq 1. This emphasizes the
point that the behavior of K as a mathematical function instead of a mathematical constant is the most fundamental and characteristic property of such mixed sample equilibria. Any polyelectrolyte properties that may be present are simply superimposed on this fundamental mixture effect. It is also important to understand the relationship of the K envelope function to the stability constants Kiof the four component ligands. It can be seen in Figure 2 that I? never becomes as large as the largest nor as small as the smallest one. This also component Ki, demonstrates the important point that K and Ki should not be used interchangeably for calculations. For natural fulvic acid samples, the difference may amount to 1or 2 orders of magnitude. Misleading AGO values for example could result from carelessness. The relative contribution of each component ligand to the K of the whole sample can be expected to change during the course of a compleximetric titration. Some important insight may be gained, therefore, from an examination of Figure 3, which shows a percent contribution curve for each of the four component ligands. The contribution curve for the most strongly binding ligand, 1, is seen to be a monotonically decreasing function of site loading. In contrast to this, the contribution curve of the weakest ligand, 4,is a monotonically increasing function of site loading. The curves for the middle ligands 2 and 3 have a very distinctive feature. They possess maximums. The positions of these maximums are indicated in the figure. Direct tests were made for the recovery of component Kivalues. The most reliable values for differential equilibrium functions can probably be obtained by careful application of computer techniques. The simple graphical method illustrated in Figure 4 serves two useful purposes, however. When quick approximate values are needed, slopes may be estimated with graph paper and a pocket calculator, as illustrated in the figure. The graphical calculation is also a simple way of demonstrating the basic nature of all differential equilibrium function calculations. At a chosen xL value, a straight line judged by eye to be tangent to the curve at that point is drawn. The calculation of a derivative always requires care, because it magnifies experimental errors. It will be seen, however, that the graphical method will give larger calculation errors than computer curve-fitting methods do. Figure 5 shows the differential equilibrium function over most of the course of the compleximetric titration. The stability constants of the four component ligands are indicated. As expected, K is found to be an envelope function. This envelope function is also expected to contain two types of K values. Original Kivalues for component ligands should be found on the curve, each with it's own xL address on the curve. The other type consists of interpolated values between the addresses of the component Kivalues. K2 and K3 have been recovered in this calculation, while K , and K4 have not. No quantitative theory is yet available for predicting this result. There are instead two empirical observations and some qualitative comments. The first observation is that the two ligands for which Kivalues were recovered were those having maximums in their percent contribution curves in Figure 3. Constants were not recovered for those lacking maximums in the percent contribution curves. The second observation is that the xL addresses of the Kivalues are at or near the maximums of the contribution curves. All four cases tested gave the same results. They are reported in Table 11. As a titration proceeds, the participation of a given ligand in the complexing has interferences from both the preceding and Environ. Sci. Technol., Vol. 22, No. 11, 1988
1331
Table 11. Recovery of Component Ki Values from the K Functions of Various Mixtures of the Four Ligands
K2 ~
mixture no.
i=1
I I1 I11 IV
1.250 0.300 0.300 1.600
scan titration point envelope no. address.xL
ligand X lo3 i=2 i=3
i =4
1.250 0.400
1.250 3.800
26 6
2.800
12
1.500
29
1.250 0.500 0.500 0.800
1.400
1.100
0.48705 0.88243 0.76303 0.43020
20.018.016.0 14.0-
12.0-
P
P x,
-
10.0
X
IY
8.06.0-
4.0
Ii
/
2.0I
I
0.10
0.20 0.30
0.40
0.50
0.60
0.70
0.80 0.90
1.00
XL I
0.90
I
I
0.80 0.70
I
I
I
I
0.60
0.50
0.40
0.30
I
I
0.20 0.10
I
0.00
xc Flgure 4. Recovery of K, values from the weighted average equilibrium function K . The envelope function K x , is used.
-
40.0
K.
-
30.0
20.0
-
Y
10.0
K4
p i I
1
0.10
0.20
0.30
I
0.40
I
I
0.50
0.60 0.70
0.80
0.90
1.00
XL
Figure 5. Differential equilibrium function, shown as an envelope function with the component constants of Table I.
following ligands. The percent contribution of the given ligand increases as the preceding ligands are increasingly used up. It decreases because of the following ligands being brought more into action. These are opposing trends. While the effects of the competing ligands cannot be expected to separately all vanish, it is reasonable to suggest 1332
Environ. Sci. Technol., Vol. 22, No. 11, 1988
K3
rec value X
‘70
lo-*
error
8.183 9.273 8.437 8.234
-1.6 +12.8 +1.4 -1.0
scan titration point envelope no. address.xL 48 16 29 46
0.16070 0.70117 0.46195 0.18753
ret
value
X
%
lo4
error
1.359 1.294 1.370 1.362
-1.5 -6.5 -0.70 -1.3
that their net collective effect might vanish through cancellation. When this happens, the maximum contribution to K and the closest approach to the Kirecovery should occur together at the point of cancellation. From this it is postulated that the prerequisite for the recovery of a component K, is that its percent recovery curve must have a maximum. If is does, then the Kican be recovered on the K envelope function at the XL location of the maximum in the contribution curve. Some important points about the addresses of Kiconstants on the K = f ( x L )envelope function should be recognized. One of the important reasons for using a known synthetic mixture is that the Ki values are known in advance. This permits the xL addresses of the Ki constants on the envelope function to be identified. The xL addresses and the Kivalues may then be used together for investigating the correctness and practical limitations of the weighted average equilibrium function theory. Once this has been done, the theory may be applied to natural mixtures with greater confidence and effectiveness. There is no available method, however, for determining the xLaddresses of components in unknown mixtures. That would only be a problem if the objective were to identify chemical species in the sample. If the sample contains a large number of components that are energetically closely spaced however, the problem becomes unimportant in any case. The reports and conclusions of other authors about the alternate methods provide some useful comparisons. According to Klotz and Hunston ( 4 ) and to Perdue (5), Scatchard type plots should not be used at all for mixtures of several components. One reason is that the nonlinear plots obtained for ligand mixtures produce only arbitrary, empirical parameters. The parameters have no chemical meaning. Saar and Weber (33) demonstrated experimentally that an isolated point value is a grossly inadequate representation of the weighted average equilibrium function. The problem is that there are no scans across the chemical composition variables. Without scans, no differential equilibrium functions can be calculated. After careful investigations of the Laplace transform method for the extraction of Kidistribution curves from K data, Shuman concluded that calculation errors seriously limit its usefulness for practical samples having large numbers of components. Theory requires the use of several consecutive differentiation steps. Theoretical distribution curves, or “affinity spectra”, are histograms of numbers of binding sites plotted against log K. With mixtures having only a few components, the calculation distorts histogram bars into smooth peaks. For samples having large numbers of components, the progressive accumulation of measurement and calculation errors can produce mathematical artifacts. This tends to confuse the chemical interpretation of the data. This is in contrast to the increased usefulness of the weighted average method for samples containing large numbers of components. The observation that the stability constants of components will be recovered on the K envelope function only at addresses defined by maximums in the percent con-
Table 111. Component Constants Recovered for Mixture IV with the Use of a Presuppressor"
scan titration point no. 19 33 49 0.354 324 0.151 420 envelope address-xL 0.621 236 rec value 3.521 X lo5 8.053 X lo4 1.367 X lo4 90 error -0.77 -3.2 -0.94 asuppressor: KaP= 1.259 X los; cL = 1.44 X 10-3. c, = 0.72 x 10-3. c, = 1.35 x 10-3.
CL
= 0.99
30 0
X
tribution curves has an important implication. It may be postulated that K1 can be recovered from an envelope function by the experimental addition of a suppressor ligand having a stability constant substantially greater than K1. By preferentially binding the Cu(I1) initially titrated in, the suppressor ligand should force the percent contribution curve of ligand 1to start a t a alow value, so that it will go through a maximum. This maximum should then define the address of K1 on the K envelope function. The ligand a-alanine has a Cu(I1) complexing stability constant of K = 1.259 X lo3, which is somewhat more than 2 orders of magnitude greater than the value of Kl. It was therefore chosen as a supressor ligand. K and K envelope functions were then calculated for mixture IV, assuming a 10 mol % addition of this supressor ligand. As expected, the numerical values of K1, K,, and K3 were all recovered from the K envelope function. The results of the calculations have been tabulated in Table 111. This result has an important implication for the practical analytical chemistry of natural fulvic acid samples. The addition of a suitable supressor ligand to a fulvic acid sample might make it possible to determine its strongest Cu(I1) stability constant. The practical significance of this is that under natural field conditions, there may well be trace level loadings of fulvic acid complexing sites by such metal ions as Cu(I1). Under those conditions, the strongest stability constants in the mixture should be manifest. The weighted average equilibrium function theory was originally developed for samples assumed to have an infinite number of components, each accounting for an infinitesimal portion of the whole sample (22). Tests with small, finite numbers of components consequently impose extreme and unrealistic strains on the theory. The most serious distortion is that the K envelope functions will have only a few widely spaced Ki values of the actual components on its curve, with most of the curve consisting of interpolation values. This makes the distinction between component Ki and interpolation K values very important. This might lead to some misunderstandings about the nature and purpose of the original theory. In particular, there may be some question about the practical usefulness of the differential K calculations. A sample with 100 ligands can give a much more favorable impression, however. If for example, the K ivalues are closely spaced and each component accounts for only 1mol % of the mixture, then the distinction between component and interpolated values becomes much less important. For an infinite number of closely spaced components, it would vanish altogether. It was necessary for these reasons to investigate the hypothetical case of a sample containing 100 component ligands, each of which accounts for 1mol % of the sample. Figure 6 shows the collection of stability constants chosen for the simulation. Ligands 1and 100 are the same as ligands 1 and 4 chosen for the four-component cases. The other values represent a hypothetical distribution. The recovery of original component Ki values from the differential K envelope function was the important test conducted with
w
b ;i 200 Y
100
1 .o
20
40
30
50
CLX103;Total ligand molarity
Flgure 6. Distribution of constants for a hypothetical mixture of ligand sites: 100 components. Vertical dashes show the range of subsequently recovered 6 values.
30.0
I
I
I
I
I
I
I
I
I
I
01
02
03
04
05
06
07
08
09
I 10
XL
Figure 7. Differential equilibrium function for the hypothetical mixture of 100 ligand sites.
this simulation. The values actually recovered lie within the range 1.95 X lo3 I CL I3.35 X lo3, as indicated in Figure 6. The K curve itself has been plotted in Figure 7. If it is assumed that the component K, values outside this range are missing because they have no maxima in the percent contribution curves, then it follows that the addition of a suitable supressor ligand should permit the recover of those a t the high end of the curve. The curve in Figure 8 confirms the predicted effect of the supressor ligand. Values in standard Gibbs energy form can be seen extending up to and including that for ligand 1. They lie between the calculation points shown. As previously mentioned, the known value of K , in the synthetic mixture was used for locating the xLaddress of the ligand 1 on the envelope function. There is another interesting observation. The stability constant of the supressor ligand has also been recovered from the plateau Environ. Sci. Technol., Vol. 22, No. 11, 1988 1333
K
d I
3.3
1
XL
Figure 8. Component K,values recovered with the use of a supressor ligand. Presented in AGO form for the hypothetical mixture of 100 ligand sites.
3.2 01
approaching XL = 1. The reason for this is that this stability constant is more than 2 orders of manitude greater than K1. If this were not so, then the value of the supressor K would be lost in just the same way that K1 is lost in the absence of a supressor ligand. This observation might be useful for the interpretation of experimental K envelope functions. Figure 8 reveals a new question that invites future investigation. The question is whether or not a supressor ligand method of analysis could be developed for simultaneously determining Kl and the xLaddress for the most reactive ligand in an unknown mixture. Figures 4 and 6 demonstrate an important practical point. When a case study is done with a realistic number and distribution of component ligands, the distinction between component and interpolation K values may become negligable compared to the experimental errors. The curve in Figure 9 tests the effect of fractionation on the recovery of component Ki values a t the high end. The top 1/4 fraction of the 100-component mixture was chosen. The first observation is that fractionation has moved the range of recovered Ki values farther up the high end. The highest recovered value is now only -20% away from K1. Secondly, it can be seen that the curve now has a slight resemblence to the shape of the curve in Figure 8. When the logarithmic nature of Figure 9 is taken into account however, it will be recognized that the effect is very weak. The apparent reason for this is that K1 differs very little from the K ivalues immediately below it. Four conclusions may be drawn from the simulations with supressor ligands. (1)The recovery of the K1 value, demonstrated in Table I11 and Figure 8, verifies the postulate that the Ki of the ith component will have an address on the K = f ( x L ) envelope function if its percent contribution curve has a maximum. (2) The recovery of the K1 value in the two tests with supressor ligand further documents the mathematical correctness of the differential K calculation. (3) A comparison of Figures 8 and 9 reveals the fact that the K1 value of the strongest complexing ligand can be recovered from an appropriate graphical plot, if it is 2 or more orders of magnitude greater than any of the other component ligands. If the difference is too small, the largest K1 cannot be recovered in this way. (4) The use of supressor ligands is potentially an important analytical chemical technique, because trace level loading of complexing sites may commonly exist under practical field conditions. The experimental K obtained will then strongly reflect the Kl of the natural ligand mixture. 1334
Environ. Sci. Technol., Vol. 22, No. 11, 1988
0.2
0.3
04
0.5
0.6
0.7
0.8
0.9
10
XL
Figure 9. Effect of fractionation on the differential equilibrium function: top ‘ I 4of the hypothetical mixture of 100 ligand sites.
The geochemical and biochemical streams of theoretical research are found to be in close agreement on basic points. Both agree, for example, that the mole fraction of unreacted binding sites defined for a whole mixture is a weighted average of the mole fractions of unreacted components of the mixture. Likewise, Klotz and Hunston have proposed the general case of eq 16 for (k),in which y is a small integer. The theory developed in the geochemical work had previously deduced the same relationship in the eq 15, but in addition, the y of the Klotz-Hunston equation is found to have a particular value y = 1. Both the Laplace transform method and the differential equilibrium method have theoretical bases for the recovery of the equilibrium constants of components from the K functions of whole mixtures. The latter is less vulnerable to the propagation of errors however, and in addition is able to account for chemical composition variables. Some doubt has been expressed in the recent literature about whether a differential equilibrium function can yield equilibrium constants that represent real chemical species. The computer simulation results have now provided a concrete demonstration that the constants of chemically real components of a mixture can be obtained. In contrast to this result, Sposito has suggested that the theoretical models found in the literature will generally produce constants that represent “quasiparticle” values (6). There has been some discussion of the meaning of equilibrium functions calculated for mixtures. The question may be resolved by examining four types of cases. (1)Integral functions calculated macroscopically for a whole unknown mixture: it is generally recognized that the numerical values do not represent real individual chemical components. (2) Differential functions calculated for a known synthetic mixture: (a) Interpolation Values. These have been found to not directly represent real chemical components. (b) Component Values. When numerical values can be recovered, they have been found to represent real chemical components. (3) Differential functions calculated for an unknown mixture with a mathematical formalism not bused on quantitative molecular level information: (a) Interpolation Values, These cannot directly represent real chemical components. (b) Component Values. Sposito (6) has correctly reported that these are not recovered. Without information about the generic chemical types of binding
sites and their actual stoichiometries, useful empirical values are recovered that represent fictional binding sites. Sposito has mentioned especially molecular level information from spectroscopic methods. (4) Differential functions calculated for an unknown mixture with a mathematical formalism based on quantitative molecular level information: (a) Interpolation Values. As before, these do not directly represent real chemical components, but may be useful for empirical characterization of a sample. (b) Component Values. Basic requirement must be met in order for Ki values of real chemical components to exist on the differential K envelope function. The mathematical formalism must have been deduced from the logic internal to the system itself. Outside hypotheses and simplifying approximations should not be used. In particular, the direct relationship between the real chemical component constants Ki and the experimental K values must have been determined. Two other points must next be considered carefully. The first is that so few investigators have heeded Sposito’s point about the significance of molecular level analytical chemical data that it is sometimes assumed that it has not been done a t all. In fact, some work has been published for samples that have been characterized (2). Various combinations of spectroscopic, potentiometric, and organic analytical data have been accumulated by a number of collaborators, This has provided the required information about the generic chemical types of binding sites and their stoichiometries. The second point has not yet been properly explained. The use of the detailed molecular level data for a natural mixture can only produce the K, values of real chemical components on the K envelope function, but cannot permit these real values to be directly identified or matched up with their chemical components. It is important to note that the inability of an investigator to identify and assign them cannot logically prove that they do not exist. They coexist on the envelope function with the interpolation values, and in limiting cases, the interpolation values may tend to vanish. Sposito has also made an important observation in the above-mentioned review. He has pointed out that ordinary complexing or binding equilibrium experiments are not sensitive to the detailed features of trace metal ion-humic substance interactions. There are two reasons why this is generally true, one mathematical and the other experimental. The mathematical reason has been demonstrated previously (25). When a measurement is made macroscopically on a whQle mixture, the numerical value obtained is a weighted average of values for all of the components of the mixture. A fundamental mathematical characteristic of the averaging process is that it damps out extreme values and attenuates sharp differences. The average never becomes as high as the highest component value nor as low as the lowest component value. In any experimental scan across the componenb, this causes a loss of important structural detail. The experimental cause is simply the ordinary measurement errors, which can further limit resolution. An important purpose of the computer simulation study was to test and demonstrate a simple mathematical method for “undoing the averaging” so that lost spectral detail can be recovered. In this way, an experiment on a mixture can become more sensitive to the detailed features of the system. Top quality analytical chemistry is required however, because otherwise the experimental errors will destroy the benefits of undoing the averaging. In a very recent publication, Buffle and Altmann (34) described the significance and limitations of the affinity
spectrum model. They pointed out that it can provide at least two kinds of insight into the properties of a system of mixed ligand sites. Considering the first of them, the affinity spectrum is a curve of the probability p(K) of the value of K plotted against K. This affinity spectrum curve will reflect changes in the chemical composition of the system. The other piece of insight relates experimental results to the inherent nature of the sample under investigation. A macroscopic measurement on a mixture produces an average value. A given analytical chemical method may be limited to a “window”that sees only part of the spectrum. It must be recognized that the numerical value of the experimental value depends on the size and location of the window. The mathematical method imposes the most important limitation. The calculations are complex and difficult and may produce mathematical artifacts. Buffle and Alt,mann (34) concluded that the calculation of a differentiala equilibrium function would be less vulnerable to the effects of experimental error. Literature Cited Perdue, E. M.; Lytle, C. R. In Aquatic and Terrestrial Humic Materials; Christman, R. F., Gjessing, E. T., Eds.; Ann Arbor Science: Ann Arbor, MI, 1983; Chapter XIV. Langford, C. H.; Gamble, D. S.; Underdown, A. W.; Lee, S. In Aquatic and Terrestrial Humic Materials; Christman, R. F., Gjessing, E. T., Eds.; Ann Arbor Science: Ann Arbor, MI, 1983; Chapter XI. Klotz, I. M. Acc. Chem. Res. 1974, 7, 162-168. Klotz, I. M.; Hunston, D. L. Arch. Biochem. Biophys. 1979, 193, 314-328. Perdue, E. M. In Humic Substances in Soil, Sediment, and Water. Geochemistry,Insolation, and Characterization; Aiken, G. R., McKnight, D. M., Wershaw, R. L., MacCarthy, P., Eds.; Wiley: New York, 1985; Chapter XX. Sposito, G. CRC Crit. Rev. Environ. Control 1986, 16, 193-229. Marinsky, J. A,; Gupta, S.; Schindler, P. J. Colloid Interface Sci. 1982, 89, 412-426. Marinsky, J. A.; Gupta, S.; Schindler, P. J. Colloid Interface Sci. 1982, 89, 401-411. Simms, H. S. J. Am. Chem. Soc. 1926,48, 1239-1261. Beckwith, R. S. Nature (London) 1959, 784, 745-746. Broadbent, F. E.; Bradford, G. R. Soil Sci. 1952,74,447-457. Schnitzer, M.; Skinner, S. I. M. Isotopes and Radiation in Soil Organic Matter Studies; International Atomic Energy Agency: Vienna, 1968; p 14. Khanna, S. S.; Stevenson, F. J. Soil Sci. 1962,93,298-305. Coleman, N. T.; McClung, A. C.; Moore, D. P. Science (Washington, D.C.) 1956, 123, 330-331. Schnitzer, M.; Desjardin, J. G. Soil Sci. Soc. Am. Proc. 1962, 26, 362-365. Schnitzer, M.; Skinner, S. I. M. Soil Sci. 1963, 96, 86-93. Himes, F. L.; Barber, S. A. Soil Sci. Soc. Am. Proc. 1957, 21, 368-373. Schnnitzer, M.; Skinner, S. I. M. Soil Sci. 1965,99,278-284. Schnitzer, M.; Hanson, E. M. Soil Sci. 1970,109,333-340. Geering, H. R.; Hodgson, J. F. Soil Sci. Soc. Am. Proc. 1969, 33,54-49. Klotz, I. M.; Hunston, D. L. J. Biol. Chem. 1984, 259, 10060-10062. Shuman, M. S.; Collins, B. J.; Fitzgerald, P. J.; Olson, D. L. In Aquatic and Terrestrial Humic Materials;Christman, R. F., Gjessing, E. T., Eds.; Ann Arbor Science: Ann Arbor, MI, 1983; Chapter XVII. Thompson, C. J.;Klotz, I. M. Arch.. Biochem. Biophys. 1971, 147, 178-185. Klotz, I. M.; Hunston, D. L. Proc. Natl. Acad. Sci. U.S.A. 1977, 74, 4959-4963. Gamble, D. S. Can. J . Chem. 1970,48, 2662-2669. Klotz, I. M.; Hunston, D. L. Biochemistry 1971, 10, 3065-3069. Cabaniss, S. F.; Shuman, S. E.; Collins, B. J. In Complexation of Trace Metals in Natural Waters; Kramer, C. J. Environ. Sci. Technol., Vol. 22,
No. 11, 1988 1335
Environ. Sci. Techno/. 1988,22,1336-1347
M., Duinker,J. C., Eds.; Martinus Nijhof: The Hague, 1984; pp 165-179.
Hunston, D. L. Anal. Biochem. 1975, 63, 99-109. Ferry, J. D. Viscoelastic Properties of Polymers, 2nd ed.; Wiley: New York, 1970. Gamble, D. S. Can. J . Chem. 1972,50, 2680-2690. Ringbom, A. Complexation in Analytical Chemistry;Interscience: New York, 1963. Martell, A. E.; Smith, R. M. Stability Constants Supple-
ment No. 1;Special Publication 25; The Chemical Society: London, 1971. (33) Saar,R. A.; Weber, J. H. Can. J.Chem. 1979,57,1264-1268. (34) Buffle, J.; Altman, R. S. In Aquatic Surface Chemistry; Stumm, W., Ed.; Wiley Interscience: New York, 1987; Chapter XIII. Received f o r review September 11, 1987. Revised manuscript received February 12, 1988. Accepted May 25, 1988.
Mathematical Modeling of the Formation of Nitrogen-Containing Pollutants. 2. Evaluation of the Effect of Emission Controls Armistead G. Russell Department of Mechanical Engineering, Carnegie Mellon University, Pittsburgh, Pennsylvania 15213
Kenneth F. McCue and Glen R. Cass" Environmental Engineering Science Department and Environmental Quality Laboratory, California Institute of Technology, Pasadena, California 91 125 W
A grid-based Eulerian airshed model is used to study
the effect of specific emission control measures on ambient NO2,total inorganic nitrate (TN), HNO,, aerosol nitrate, PAN, NH,, and ozone concentrations in the Los Angeles area. NO, and reactive hydrocarbon (RHC) emission reductions of up to 61 and 37%, respectively, are examined. NO2 and T N concentration reductions in excess of 50% averaged over 20 monitoring sites are achieved at the highest level of emission control studied. The distribution of T N air quality improvements between HNO, and aerosol nitrate is affected by the NH, emission rate of the NO, control technologies employed. Peak 1-h O3 concentrations at many sites in the eastern portion of the air basin studied decline by more than 25% at the highest NO, and RHC control levels studied. Introduction In part 1 of this series, the performance of a grid-based
photochemical airshed model for NO2, total inorganic nitrate (TN), PAN, HNO,, NH,, aerosol nitrate (AN), and ozone formation and transport was evaluated (I). Model predictions were compared against experimental observations made for this purpose in the Los Angeles area over the period August 30-31,1982. It was found that O3 and PAN concentration predictions were in excellent agreement with observations and that NO2 predictions were in closer agreement with observed values than in many previous studies. On average, TN, NH,, and HNO, concentration predictions differed from observations by very small absolute amounts: 2.7 pg rn-, (1.1ppb), 0.7 ppb; and 4.2 pg rn-, (1.65 ppb), respectively. In this paper, that model will be used to explore the effect of specific emission control measures on the HN03/NH3/aerosol nitrate system. Newly adopted ambient air quality standards for particulate matter in sizes less than 10-pm diameter (PM,,) have focused attention on the importance of aerosol nitrate control. Recent measurements of PMlo concentrations and chemical composition in the Los Angeles area show a peak 24-h average PMlo concentration of 299 pg/m3 at Rubidoux near Riverside, CA, during the year 1986 (2). That PMlo concentration is approximately twice as high as allowed by the newly adopted PMlo standards of the U.S. Environmental Protection Agency, and 6 times higher than the State of California 24-h average PMlo ambient air quality standard. On the peak day in question (October 29, 19861, aerosol 1336
Environ. Sci. Technol., Vol. 22, No. 1 1 , 1988
nitrate plus associated ammonium ion accounted for more than 40% of the observed PMlo peak value (2). Clearly, it will not be possible to attain compliance with PMlo air quality standards in the Riverside area without considering the approaches available to reduce aerosol nitrate concentrations via emission controls. Key questions that must be addressed include: (1)What effects would controls on reactive hydrocarbons (RHC) and oxides of nitrogen (NO,) emissions have on aerosol nitrate concentrations? (2) How would controls on ammonia emissions affect HNO, and NH4N03concentrations? (3) If NH, injection technology is used for NO, control at stationary sources, would the NH, emissions from such systems aggravate the ",NO, aerosol control problem? While the principal focus of this work is directed a t examining the effect of emission controls on the HNO3/NH3/NH,NO3 system, further information on the response of NOz, O,, and PAN concentrations to emission controls is obtained as a byproduct of the analysis. Emission Control Opportunities
Emission control measures evaluated as part of this study are itemized in Table I. That table has been divided into five groups. Group 1controls reflect a subset of the reduction possibilities that have been documented as part of the 1982 Air Quality Management Plan (AQMP) for the South Coast Air Basin that surrounds Los Angeles (3). This group of controls approximates the effect of many of the emission reductions that were expected to be implemented in the Los Angeles area in the years following the 1982 base year, but without extension of vehicular catalyst utilization or ammonia injection technology beyond that used in 1982. Group 2 and group 3 controls simulate the effect of fleetwide improvements in emissions from motor vehicles at target levels that have been discussed by state and federal regulatory agencies (7, 9). Group 4 and group 5 controls would further reduce NO, emissions from stationary sources through the use of noncatalytic ammonia injection or selective catalytic reduction (SCR) technology. The 1982 emission inventory employed during the model verification effort of part 1 of this study ( I , 15) will be referred to as the base case. The 1982 base case emissions from each source class that will be considered for control are given in Table I, along with the percentage reduction in those emissions that would result if the control measures
0013-936X/88/0922-1336$01.50/0
0 1988 American Chemical Society
