Research: Science and Education
How To Compute Labile Metal–Ligand Equilibria Robert de Levie Department of Chemistry, Bowdoin College, Brunswick, ME 04011;
[email protected] Insofar as metal–ligand complexes are substitution-labile, the formal description of their speciation resembles that of acid–base problems, which usually involve sufficiently fast kinetics that, on the time scale of seconds, can be described in terms of chemical equilibria. However, with the exception of chelating ligands that predominantly form one-to-one complexes with metal ions, treatment of these equilibria is gradually disappearing from undergraduate textbooks, in order to make room for topics that take less space or are considered to be more relevant. The relevancy issue will not be considered here, since this is largely a matter of taste, which falls outside the realm of meaningful discussion. Here we will, instead, consider the matter of space taken up by the traditional treatment of such equilibria, in terms of successive approximations. The chemical complexity of many complexation equilibria is, of course, the main reason for their seeming intractability. Traditionally, students either are given a black-box approach using commercial, “canned” programs or are presented with a set of mass balance and charge balance equations, plus the expressions for the relevant equilibrium constants. In the latter case they are then guided to an answer by using a trial-and-error approach. Both methods can of course yield correct answers. In didactic terms, however, use of a canned program is non-transparent. Moreover, because it involves specialized hardware, is not readily transferable. On the other hand, the trial-and-error method is uninspiring and time-consuming in class and space-consuming in textbooks. Yet, chemical complexity per se is not an academic exercise in futility but a fact of life, and seeing how simple chemical relations for individual steps indeed apply in complicated systems illustrates an important idea underlying science, viz., that the over-all complexity we see all around us can be understood in terms of the interplay between many much simpler steps. Hence the question: Can the trial-and-error treatment be streamlined, while remaining transparent to the user? The answer is “yes”, and this article will indicate how. One can simplify the approach by making the computer do its tedious, numerical part. Many students may already be familiar with the iterative method used in, for example, nonlinear least squares, which they may have encountered on their laptop or desktop computer in the form of Excel’s Solver, a generalized reduced gradient alogrithm developed by Lansdon and Waren (1) and implemented in Excel by Fylstra and Watson (2). It has counterparts in virtually all numerical software packages and is typically associated with the names of Levenberg (3) and Marquardt (4). Students with some exposure to numerical analysis may know its single-parameter predecessors, such as bisection or the Newton– Raphson method, the latter still surviving in Excel as the GoalSeek method, or may be familiar with the simplex method of Nelder and Mead (5) or any of the other optimization algorithms discussed, for example, in chapters 9 and 10 of the Numerical Recipes (6). All of these are examples of 136
Journal of Chemical Education
•
iterative methods for finding numerical answers to mathematically well-posed problems, methods that are readily available on many computers, and that many students already know. The task of the teacher, then, is to formulate the problem in such a way that it becomes suitable for such an iterative computer solution. The purpose of this article is to illustrate how this can be done. The Ligand Function Our point of departure will be the ligand function LF, constructed by analogy with the proton function HF (7) to represent the concentrations of “gainers” and “losers”, respectively, of ligand, as measured with respect to the original species added. For example, when one dissolves CdCl2 in water at a concentration of C1 moles per liter (acidified with, for example, C2 M HNO3 in order to avoid hydrolysis), and recognizes that Cd(II) can form mononuclear chloro complexes up to CdCl42−, one can write Cl
Cl −
F =
+ 2 CdCl 4 2 − + CdCl 3 − − CdCl + − 2 Cd 2 +
= 0
(1)
where water has “gained” chloride as Cl−, CdCl42− has gained two Cl− ligands, CdCl3− has gained one Cl−, CdCl+ has lost one Cl−, and Cd2− has lost two Cl−, all versus CdCl2. This result can be written simply by inspection, as just indicated. Alternatively, eq 1 can be derived from the relevant mass and charge balance expressions, which in this case read
Cd 2 +
+ CdCl +
+ CdCl 2 + CdCl 4 2 −
+ CdCl 3 − Cl −
+ CdCl +
+ 4 CdCl 4 2 −
− OH −
NO3 − 2 Cd 2 +
(2)
+ 2 CdCl 2
+ 3 CdCl 3 − H+
= C1
+ CdCl +
= 2C1
(4)
= C2
(5)
= C2
+ H+
− CdCl 3 −
(3)
− Cl − − 2 CdCl 4 2 −
− NO3 − − OH −
(6) = 0
The chloride function, eq 1, focuses narrowly on chloride, to the exclusion of spectator ions such as [H+], [NO3−], and [OH−] that are peripheral to the problem posed.
Vol. 84 No. 1 January 2007
•
www.JCE.DivCHED.org
Research: Science and Education
The numerical solution of eq 1 is readily obtained as follows. First one rewrites it in terms of the total analytical cadmium chloride concentration C1 and the associated equilibrium constants as Cl
F =
Cl −
+ 2 CdCl 4 2 −
+ CdCl 3 −
− CdCl + =
Cl −
+ 2α 4 + α 3 − α1 − 2α 0 C1 2β 4 Cl −
=
Cl −
− 2 Cd 2 +
4
+ β3 Cl −
− β1 Cl +
β 4 Cl −
4
−
+ β3 Cl −
+ β 2 Cl −
2
3
(7)
− 2 C1 = 0
3
+ β1 Cl − + 1
using the traditional formalism of concentration fractions α and overall formation constants β for metal–ligand complexes, that is,
CdCl +
= α1C1 = β1 Cd 2 + Cl −
CdCl 2 = α 2C1 = β 2 Cd 2 +
Cl −
(8) 2
CdCl 3 −
= α 3C1 = β3 Cd 2 +
Cl −
CdCl 4 2 −
= α 4C1 = β 4 Cd 2 +
Cl −
(9)
3
4
(10) (11)
After substitution of the numerical values for these equilibrium constants, one then lets Excel’s Solver (or GoalSeek, or their equivalents in other software, or even a simple bisection routine) adjust [Cl−] so that ClF becomes zero. Once [Cl−] is known, the concentrations of all metal–ligand complexes follow directly from eq 2 plus expressions such as eqs 8 through 11. All this is illustrated with Excel in Figure 1, and its results can be compared directly with those given in ref 8. As can be seen, once the problem has been set up by entering the equilibrium constants in cells B1:B4, the total analytical concentration C1 in B5, an initial estimate for [Cl᎑] in cell B7 and defining the chloride function ClF in cell B8, the numerical solution simply consists of finding the root of a single equation, once. Cells E1:E5 then automatically compute the concentrations of the various cadmium species, and cells E7:E9 verify that the mass and charge balance equations are indeed correct. The computed results are shown in Figure 1 to many more decimal places than justified in terms of the limited experimental precision of the equilibrium constants. A higher numerical precision is displayed here only to illustrate more clearly the principle of the method. The method per se is not approximate, and its numerical answer will be as precise as the software allows. If so desired, one can then use the spreadsheet to construct the corresponding logarithmic concentration diagrams, or other aids to visualize such a complex equilibrium. The formulas for doing so can be found in Figure 1, where they www.JCE.DivCHED.org
•
Figure 1. The calculation of the speciation of cadmium chlorides in acidified CdCl2 (total analytical concentration 0.01 M) on an Excel spreadsheet. The overall formation constants (taken from ref 8) as log β1 = 1.98, log β2 = 2.6, log β3 = 2.4, and log β4 = 1.7, are stored in cells B1:B4, and the chosen value of C1 in B5. Cell B7 initially contained an estimated value, and was refined by Solver to the value shown, after cell B8 (containing ClF) had been designated as the target cell, to be minimized by changing the value of [Cl−] in cell B7. The data in E1:E5 compute the concentrations of the cadmium species, while E7 and E8 test the mass balances of cadmium and chloride, respectively, and E9 the charge balance. The data are shown with many more decimal places than justified by the experimental precision of the input data, merely to show the achievable numerical precision.
Figure 2. The chloride function ClF vs pCl for the parameters shown in Figure 1. The location of the root of this function (where it crosses zero) is emphasized by a gray-filled circle. This root can be found, for example, with Excel’s Solver, GoalSeek, or by bisection.
Vol. 84 No. 1 January 2007
•
Journal of Chemical Education
137
Research: Science and Education
are listed specifically for cells E1:E5, that is, for the concentrations of the various cadmium species, in terms of the chloride concentration in cell B7. How easy is it for a computer to find this solution? Figure 2 illustrates that the ligand function LF is a continuous, monotonically decreasing function of pL, just as the proton function HF is of pH. As can be seen from eq 7, for [Cl−] → +∞, that is, pCl → ᎑∞, one finds that ClF → [Cl−] → +∞. On the other hand, for [Cl⫺] → 0, pCl → ᎑∞, and one has Cl F → ᎑2[CdCl42⫺] → ᎑2C1. In our example, C1 = 0.01, hence the lower limit of ClF is ᎑0.02 M. The root (i.e., zero-crossing) of ClF can therefore be determined readily with a number of single-parameter or multiple-parameter root-finding methods, available in almost any numerical software package. We use Excel’s Solver here merely because of its wide availability.
as an explicit function of [Cl−] by Vt = − Vs
d LF d LF = d [L ] d ln [L ]
(12)
which can be used to describe the extent to which the solution resists a change in ligand concentration [L] and which can be compared directly with the proton buffer strength H
B = [ H]
d HF d HF = d ln [H] d [ H]
(13)
The same formalism also provides simple expressions for the titration of metal ions with their ligands. In this case we have, again in complete analogy with the acid–base case,
Vs LFs = −Vt LFt
(14)
where the sample and titrant are denoted by the subscripts s and t, respectively, while the original sample volume is Vs, whereas Vt is the volume of titrant added at any particular stage during the titration. The chloride function LFs has the form of that of the original sample, and LFt that of the titrant, but the individual concentrations in these ligand functions change during the titration. This is most readily illustrated with the following, worked-out example. For the sample we have, starting from Cd(NO3)2, Cl
Fs =
Cl −
+ 4 CdCl 4 2 −
+ 3 CdCl3 −
+ 2 CdCl 2 + CdCl +
(15)
which originally is of course zero, and for the titrant Ct M NaCl we have Cl
Ft =
Cl −
(16)
− Ct
so that the entire progress of the titration can be described 138
Ft
=
=
Journal of Chemical Education
•
+ 3 CdCl 3 −
+ 2 CdCl 2 + CdCl −
+ Cl −
C t − Cl −
( 4α 4
+ 3α 3 + 2α 2 + α1 )C s + Cl − C t − Cl −
4β 4 Cl −
In the case of acid–base equilibria, the proton function is most useful in simplifying the descriptions of buffer action and of titrations. The same applies to metal–ligand equilibria. Since this communication focuses on simple equilibrium calculations, we will only briefly indicate these two aspects here. One can define a ligand buffer strength LB as
B = [L ]
Fs
Cl
4 CdCl 4 2 −
Buffer Action and Titrations
L
Cl
4
+ 3β3 Cl − 2
+ 2 β2 Cl − β 4 Cl
− 4
+ β 2 Cl −
+ β1 Cl −
+ β 3 Cl 2
3
− 3
+ β1 Cl −
=
C s + Cl −
+1
(17)
C t − Cl −
where Vt and [Cl−] vary during the titration. Equation 17 yields the ratio Vt兾Vs directly for any given value(s) of [Cl−]. Alternatively, eq 17 can be solved iteratively for [Cl−] for any given value(s) of Vt兾Vs. Incidentally, the equilibrium constants for the many metal–ligand complexes that one finds listed in compilations such as ref 9 were often obtained by fitting experimental titration data to expressions equivalent to eq 17. The general eq 14 in terms of the chloride functions leads directly to the complete expression for the progress curve of this ligand titration, just as the corresponding proton function does for an acid–base titration. The mathematical complexity of eq 17 reflects the participation of a number of species in this titration, just as it would if we were to titrate a tetraprotic base with strong monoprotic acid. The latter example is, indeed, quite analogous, except for the different conventions for equilibrium constants used in the two subdisciplines: formation constants for metal–ligand complexation, dissociation constants for acid–base equilibria. Ligand Competition The advantages of using ligand functions are most evident in chemically fairly complicated cases, as when different ligands compete for the same metal ion. This situation occurs often, especially in basic solutions, where hydroxide is a common second ligand. Here we will illustrate how such cases can be solved easily using ligand functions, again because these allow us to replace repeated, manual trial-anderror efforts by a single computer-controlled iteration procedure. As our example we will use the formation of mixed complexes of Hg2+ with both Cl− and OH−. Specifically, we will calculate the various mercury(II) species in an unbuffered aqueous solution containing a total analytical HgCl2 concentration C of 0.27 M.
Vol. 84 No. 1 January 2007
•
www.JCE.DivCHED.org
Research: Science and Education
We will make the following assumptions: (i) that Hg(II) only forms the mono-nuclear complexes listed in ref 10, that is, HgCl+ through HgCl42− plus Hg(OH)+ and Hg(OH)2 and (ii) that the water used is not contaminated (as it often is) by dissolved CO2. We now have two ligands, and therefore two ligand functions, ClF and OHF, respectively, which at equilibrium should both be zero, simultaneously: Cl
+ 2 HgCl 4 2 −
Cl −
F =
F =
OH
−
OH − 2β 2 ′ OH −
+
β 4 Cl −
(18)
+ β1′ OH −
+ β3 Cl −
3
+ β 2 Cl −
C 2
+ β1 Cl − + β 2 ′ OH −
− H+
−
= 0
+ 2 Hg (OH)2 + HgOH +
4
2
+ β1′ OH − + 1
− 2 HgOH +
− 2 Hg (OH)2 OH
F =
+ HgCl 3 −
− 2 Hg 2 +
− HgCl +
OH
= 0
(19)
2
(25)
Kw OH −
where we have used the same formalism for the mercury chlorides as used earlier for the cadmium(II) chlorides, with an apostrophe ⬘ identifying the corresponding hydroxo complexes, and where Kw = [H+][OH−] defines the autoprotolysis constant of water.
Both expressions can be written down by inspection, for example, simply by considering which species have gained chloride or hydroxide ligands, respectively, and how many. Alternatively, they can be derived from the mass balance expressions for Hg, Cl, and OH,
Hg 2 +
+ HgCl +
+ HgCl 4 2 − Cl −
+ HgCl 2 + HgCl 3 −
+ HgOH +
+ HgCl +
+ Hg (OH)2
= C
(20)
+ 2 HgCl 2 + 4 HgCl 4 2 −
+ 3 HgCl 3 −
= 2C
OH − + 2 Hg (OH)2 + HgOH + − H +
(21)
= 0 (22)
and the charge balance (which here, as in eq 6, is redundant). 2 Hg 2 +
+ + HgCl + + HgOH + H + − Cl − − HgCl 3 −
− 2 HgCl 4
2−
− OH
−
(23)
= 0
We now write eqs 18 and 19 explicitly in terms of [Cl−] and [OH−] as 2β 4 Cl −
4
+ β3 Cl −
− β1 Cl −
− 2
− 2 β2 ′ OH − Cl
3
2
− 2β1′ OH −
C
−
F = Cl +
β 4 Cl −
4
+ β 2 Cl
+ β3 Cl − − 2
+ β 2 ′ OH −
3
+ β1 Cl 2
(24) −
+ β1′ OH −
www.JCE.DivCHED.org
•
+ 1
Figure 3. The calculation of the speciation of the chloro- and hydroxo-complexes of an unbuffered aqueous solution made to contain a total analytical concentration of 0.27 M HgCl2 on an Excel spreadsheet. The overall formation constants (listed in ref 8) were used, that is, log β1 = 6.74, log β2 = 12.22, log β3 = 14.07, log β4 = 15.07, log β1⬘ = 10.30, and log β2⬘ = 21.70. In the optimization, the value of ClF 2 + OHF 2 in the target cell B14 was minimized by Solver by changing the values of [Cl⫺] and [OH⫺] in cells B10:B11.
Vol. 84 No. 1 January 2007
•
Journal of Chemical Education
139
Research: Science and Education A
B
C
Figure 4. The ligand functions for chloride and hydroxide for the example and parameter values used in Figure 3: (A) plotted as a function of pCl at constant pOH = 9.888 and (B) plotted as a function of pOH at constant pCl = 3.7234. In A and B the zero-crossings of ClF and OHF coincide and are indicated with a gray-filled circle. (C) As long as the computation has not homed in on the equilibrium condition, the ligand functions for chloride and hydroxide (with pCl kept constant at 3.0 rather than at 3.7234) have different zero crossings, as indicated by a gray circle for Cl and a black one for OH.
140
Journal of Chemical Education
•
At equilibrium, both of these ligand functions must be zero, simultaneously. We can achieve that by setting to zero the sum of their squares, that is, by letting Solver adjust both [Cl−] and [OH−] in order to set the quantity {(ClF )2 + (OHF )2 } equal to zero. (For this multiparameter adjustment one should of course use Solver or some other multiparameter minimum-finding routine.) Doing so yields both ligand concentrations in one single operation, whereupon one can calculate the concentrations of all mercury species involved. This approach is illustrated in Figure 3, and the results can be compared with those on pp 296–300 of ref 10. In order to have easily checked, round numbers for the mass and charge balances, and to stay clear of precipitation, we have used C = 0.27 M rather than [HgCl2] = 0.27 M as in ref 10. Consequently, in comparing the above results with those of ref 10, all our concentrations should be increased by 0.324%. The quantity {(ClF )2 + (OHF )2 } used here is not unique, and other parameter combinations that make both ClF and OH F approach zero simultaneously, such as {|ClF | + |OHF |}, can be used instead. Of course, when different functions weigh ClF and OHF differently, the numerical pathways to their roots will also be slightly different. A note on Solver may be helpful here. Solver is an implementation, by Frontline Systems, of the Lansdon–Waren alogrithm and is incorporated in Microsoft Excel as well as in IBM Lotus 1-2-3 and Corel Quattro Pro. For scientists it is one of the most useful spreadsheet routines. In Excel, it must be activated specifically, preferably when the program is installed. If that installation used the minimum version of Excel, one may need to go back to the original installation disk, or (for institution-wide software) ask the program administrator to do so. Solver usually handles multiple adjustable parameters adroitly, especially when using its automatic scaling to minimize problems when those adjustable parameters have quite different magnitudes. In the current example, Solver is not quite so deft even though, given the experimental uncertainties in the values of the equilibrium constants used, Solver in its standard (double-precision) implementation is quite adequate for the above task. However, to illustrate that our approach is indeed numerically correct, the results shown in Figure 3 were obtained by further refinement using a Levenberg–Marquardt alogrithm. Figures 4A and B illustrate the two ligand functions involved here, each as computed by keeping the other at its minimum value. Again, the ligand functions LF are continuous, monotonically decreasing functions of the appropriate pL when the other ligand concentration (but clearly not its ligand function) is kept constant. This latter condition pertains in routines such as Solver, which adjusts one parameter at a time. As long as the computation has not yet reached equilibrium, the two roots do not yet coincide, as illustrated in Figure 4C. The compilation of Smith and Martell (9) also lists the mixed complex HgOHCl, for which it gives three roughly equivalent formation constants, leading to β⬙ = [HgOHCl]兾([Hg2+][OH−][Cl−]) with numerical values for log β⬙ between 17.5 and 17.9. To incorporate this species, add a term ᎑[HgOHCl] to ClF, and a term +[HgOHCl] to OH F. The list of compounds in D1:E7 of Figure 3 must then be extended to include [HgOHCl], the concentration of
Vol. 84 No. 1 January 2007
•
www.JCE.DivCHED.org
Research: Science and Education
which is most readily computed as [HgOHCl] = β⬙[Hg2+][OH−][Cl−] = {β⬙C[OH−][Cl−]}兾denom, where denom is the common denominator in eqs 23 and 24 as well as in the formulas given for E1:E7 in Figure 3. This denominator must, likewise, be extended to include a term +β⬙[OH−][Cl−]. Similarly one can include polynuclear complexes such as Hg2OH3+ or Hg3(OH)33+. The inclusion of mixed or polynuclear complexes of course involves additional terms, but neither interferes with nor limits the above approach, which retains its simplicity of only needing one computer-performed optimization step, no matter how many different species are involved. Buffer Action and Titrations In situations with more than one ligand, one must also define more than one ligand buffer strength, one for each ligand. The appropriate definition in that case is a simple extension of eq 12, viz., L
L
Bi = Li
L
∂ F i ∂ L i
= j ≠ i
∂ F i ∂ ln L i
(26) j ≠ i
in terms of partial derivatives. For example, in the case of the mixed chloro and hydroxo complexes of mercury(II), one must consider two separate buffer strengths, ClB and either OH B or HB. Ligand titrations will almost always involve a deliberate change in the concentration [Li] of only one of the ligands, while buffering the other ligand concentrations [Lj]j≠i. Such a titration is completely described by the corresponding progress equation Vs LiFi , s = −Vt LiFi , t
(27)
and its routine application can often be simplified by using conditional equilibrium constants, appropriate for the specific conditions under which the titration is performed. Discussion In the treatment of labile metal–ligand complexes, the ligand function LF is the analogue of the proton function HF in acid–base behavior (7). Its calculation is readily performed by computer, making it a convenient tool to simplify the numerical evaluation of metal–ligand equilibria without any loss of rigor. The ligand function allows us to relegate otherwise tedious iterations to a computer, while retaining complete control over what is calculated.
www.JCE.DivCHED.org
•
The approach is illustrated first with an example involving only a single ligand and subsequently with a case exhibiting two competitive ligand equilibria. When a single ligand is involved, the mathematical problem can be solved either with Solver, as illustrated here, or by letting the computer find the root of the corresponding ligand function with, for example, bisection or a Newton–Raphson iteration. When more than one ligand is involved, one instead must use an appropriate combination of ligand functions that forces them simultaneously to zero. Unfortunately, neither the sum of squares nor the sum of the absolute values provides a combined function with a zero crossing, so that one must use a multiparameter method that finds a minimum. Ligand buffering and metal–ligand titrations can also be treated conveniently by using the corresponding ligand function formalism, but these are not discussed here in any detail. Acknowledgements I thank Jeffrey Nagle and Keith Oldham for their helpful comments, and especially the latter for his suggestion to minimize the sum of the absolute values of the ligand functions as an alternative to minimizing the sum of their squares. Literature Cited 1. Lasdon, L. S.; Waren, A. D. In Design and Implementation of Optimization Software; Greenberg, H. J., Ed.; Sijthoff and Noordhoff: Alpen a/d Rijn, 1978; pp 363–397. 2. Fylstra, D.; Lasdon, L.; Watson, J.; Waren, A. Interfaces 1998, 28 (5), 29–55. 3. Levenberg, K. Q. Appl. Math. 1944, 2, 164. 4. Marquardt, D. W. J. Soc. Ind. Appl. Math. 1963, 11, 431– 441. 5. Nelder, J. A.; Mead, R. Comp. J. 1965, 7, 308–313. 6. Press, W. H.; Flannery, B. P.; Teukolsky, S. A.; Vetterling, W. T. Numerical Recipes: The Art of Scientific Computing; Cambridge University Press: Cambridge, 1986. 7. de Levie, R.; Chem. Educator 2001, 6, 272–276; in this article we used HC instead of HF in order to emphasize that the proton function is a sum of concentrations. 8. Butler, J. N.; Cogley, D. R. Ionic Equilibrium: Solubility and pH Calculations; Wiley: New York, 1998; pp 245–248. 9. Smith, R. M.; Martell, A. E.; Critical Stability Constants; Plenum Press: New York, 1976; Vol. 4. 10. Butler, J. N. Ionic Equilibrium: A Mathematical Approach; Addison–Wesley, Reading, MA, 1964. 11. Software Download Home Page. http://digilander.libero.it/foxes/ SoftwareDownload.htm; see also Exellaneous Home Page. http:// www.bowdoin.edu/~rdelevie/excellaneous (both accessed Sep 2006).
Vol. 84 No. 1 January 2007
•
Journal of Chemical Education
141