An Expansion for Quasi-Reversible Linear Potential Sweep

obtained in the form of an exponential (double) series. Euler's transformation of series is used for practical evaluation. Numerical precision problem...
4 downloads 0 Views 136KB Size
Correspondence Anal. Chem. 1997, 69, 1465-1469

An Expansion for Quasi-Reversible Linear Potential Sweep Voltammetry and the Use of Euler’s Transformation of Series Davis K. Cope

Department of Mathematics, North Dakota State University, Fargo, North Dakota 58105-5075

An expansion for the current-voltage curves occurring in quasi-reversible linear potential sweep voltammetry is obtained in the form of an exponential (double) series. Euler’s transformation of series is used for practical evaluation. Numerical precision problems limiting past applications of Euler’s transformation can be overcome by the use of arbitrary precision arithmetic, now commonly available in mathematical software. Sample calculations are given, and the relation of the solution to previous work in the area is described. This paper reports a new series solution for quasi-reversible linear potential sweep voltammetry and its practical evaluation by Euler’s transformation of series. Randles1 and Sevcˇik2 obtained the first theoretical results for linear potential sweep voltammetry by determining the currentvoltage curves for a reversible electrode process. In this case, the curves can be reduced to a single dimensionless curve, calculated by Randles using a numerical/graphical method and by Sevcˇik from an integral representation. Delahay3 formulated linear potential sweep for an irreversible electrode process. The current-voltage curves can also be reduced to a single dimensionless curve, which Delahay obtained by numerical solution of an integral equation. Matsuda and Ayabe4 formulated linear potential sweep for quasi-reversible electrode processes and obtained the current-voltage curves by numerically solving an integral equation. The quasi-reversible case differs from the reversible and irreversible cases in that its solution does not reduce to a single curve. Even after maximum simplification by the introduction of dimensionless quantities, the full solution requires a two-parameter family of curves with parameters R (the transfer coefficient) and Λ (a dimensionless sweep rate; see eq 1). Methods for the generation of such curves are therefore of particular interest. Reinmuth5 found new representations in terms of exponential series expansions for the solutions of the reversible E-mail: [email protected]. (1) Randles, J. E. B. Trans. Faraday Soc. 1948, 44, 327-338. (2) Sevcˇik, A. Collect. Czech. Chem. Commun. 1948, 13, 349-377. (3) Delahay, P. J. Am. Chem. Soc. 1953, 75, 1190-1196. (4) (a) Matsuda, H.; Ayabe, Y. Z. Elektrochem. 1955, 59, 494-503. (b) Bard, A. J.; Faulkner, L. Electrochemical Methods; Wiley: New York, 1980; Chapter 6. (5) (a) Reinmuth, W. H. Anal. Chem. 1961, 33, 1793-1794. (b) Reinmuth, W. H. Anal. Chem. 1962, 34, 1446-1454. S0003-2700(96)00826-8 CCC: $14.00

© 1997 American Chemical Society

and irreversible cases. Further work on linear potential sweep voltammetry includes a remarkably efficient series solution for the reversible case obtained by Oldham6 using a semi-integration approach and an exponential series solution for the quasireversible case obtained by Myland and Oldham7 under the restriction that the transfer coefficient R is rational. In the first part of this paper, we begin with a series expansion in reciprocal powers of Λ given by Matsuda and Ayabe.4 In its original form, this expansion is too complicated for actual calculations, and Matsuda and Ayabe only used it as a means of unifying discussion of the reversible (Λ ) ∞), irreversible (Λ ) 0), and quasi-reversible (0 < Λ < ∞) cases. Here we show that this expansion can be developed into an exponential (double) series, a form usable for practical calculations. Consequently, exponential series expansions are now known for all basic electrode processes: reversible, irreversible, and quasi-reversible. The present expansion also appears closely related to the expansion of Myland and Oldham,7 but the restriction to rational R is removed. Reinmuth’s exponential series for the reversible case only converges for a portion of the potential sweep range, and the same difficulty occurs with the expansion obtained here. The use of Euler’s transformation of series to extend the range of convergence was suggested by Reinmuth5 but did not prove successful: “attempts to use that approach in this laboratory ... produced erratic results”.8 In the second part of this paper, a numerical instability problem associated with the calculation of coefficients for Euler’s transformation is described. This problem is the likely source for previous difficulties in implementing Euler’s transformation, and it can be substantially alleviated by using the arbitrary precision arithmetic now available in many types of computational software. We apply Euler’s transformation to the present case, use arbitrary precision arithmetic for calculation of the transformation, find that the instability problem is adequately overcome and that the range of convergence is indeed extended, and present the corresponding numerical results for quasi-reversible linear potential sweep voltammetry. (6) (a) Oldham, K. B. J. Electroanal. Chem. 1979, 105, 373-375. (b) Oldham, K. B. SIAM J. Math. Anal. 1983, 14, 974-980. (7) Myland, J. C.; Oldham, K. B. Anal. Chem. 1994, 66, 1866-1872. (8) Nicholson, R. S.; Shain, I. Anal. Chem. 1964, 36, 706-723.

Analytical Chemistry, Vol. 69, No. 7, April 1, 1997 1465

Since the reversible reaction expansion is a special case of the present quasi-reversible reaction expansion, these results validate Reinmuth’s original suggestion5 regarding the use of Euler’s transformation.



f(u) du

x

xπ(x - u)

) g(x)

-∞

has the solution DESCRIPTION In linear potential sweep voltammetry, soluble species O and R in dilute concentration in a solution containing excess supporting electrolyte diffuse to a shielded planar electrode with area A and undergo simple electron transfer O + nee- h R in response to an applied potential E ) E0 - υt, where υ is sweep rate (>0) and t is time. The diffusion coefficients for O and R are DO and DR; the surface concentrations at the electrode are (CO)s, (CR)s; the bulk concentration for O is C*O; it is assumed that R is initially absent and has zero bulk concentration. Designating the reversible and quasi-reversible currents by IRV and IQR, respectively, the reversible case is determined by the Nernst relation,

f(x) )

d dx



g(u) du

x

xπ(x - u)

-∞

Indeed, eq 2a can be solved immediately:

iRV(x) )



d dx

x

-∞

exp(u) du exp(u) + 1 xπ(x - u)



)

∑ (-1)

xm + 1 exp((m + 1)x)

m

(3a)

(3b)

m)0

(CR)s/(CO)s ) exp(ne f (E°′ - E))

and the quasi-reversible case by the Butler-Volmer relation,

IQR/ne FA ) k°′ exp[-Rne f (E - E°′)](CO)s k°′ exp[(1 - R)ne f (E - E°′)](CR)s

where E°′ is the formal potential, k°′ is the corresponding formal heterogeneous rate constant, R is the transfer coefficient, and f ) F/RT has its usual meaning. Introducing the dimensionless quantities

The integral solution is equivalent to that obtained by Sevcˇik,2 and the series is the solution obtained by Reinmuth.5 Notice the series only converges for x < 0. By the same procedure, an alternate form for the quasi-reversible eq 2b can be obtained:

iQR(x) +

1 d Λ dx



x

-∞

exp((1 - R)u) iQR(u) du ) iRV(x) (4) exp(u) + 1 xπ(x - u)

GENERAL EXPANSION Matsuda and Ayabe4 obtained the following expansion for the quasi-reversible current:

( )



DR 1 ln x ) ne f (E1/2 - E) with E1/2 ) E°′ + 2ne f DO

iRV(u) du

x

xπ(x - u)



iQR(u) du

x

xπ(x - u)

-∞

+

)

exp(x) exp(x) + 1

(1)

(2a)

exp(x) 1 exp((1 - R)x) iQR(x) ) Λ exp(x) + 1 exp(x) + 1 (2b)

The lower bound x ) -∞ is a mathematical convenience corresponding to an idealized experiment; the resulting solutions show negligible current until the potential enters a significant range. In practice, of course, the experiment starts at a definite time at a potential of no detectible current and is swept through a physically significant range. Both integral equations are closely related to the Abelian integral equation, which has a well-known explicit solution: 1466

Analytical Chemistry, Vol. 69, No. 7, April 1, 1997

in(x)

i0(x) ) iRV(x)

it can be shown4,7 that the dimensionless currents iRV and iQR satisfy the integral equations

-∞

-n

(5)

where

iQR(x) ) IQR/(neFAC*O(ne f υDO)1/2)



n

n)0

1/2 iRV(x) ) IRV/(neFAC* ) O(ne f υDO)

R 1/2 Λ ) k°′/(ne f υD1-R O DR)

∑(-1) Λ

iQR(x) )

in(x) )

d dx



x

-∞

(6a)

exp((1 - R)u) du (6b) i (u) exp(u) + 1 n-1 xπ(x - u)

This result is obtained by substituting expansion 5 into the integral eq 4 and matching powers of Λ. It should be observed that i0(x) is independent of R and, for n > 0, that each in(x) depends on R but not Λ. Consequently, as Λ f ∞, series 5 for the quasireversible current iQR(x) reduces to the single term i0(x), that is, to the reversible current iRV(x), as one would expect. Matsuda and Ayabe did not use expansion 5 for computation, turning instead to numerical methods for the integral eq 2b. Here we obtain the following exponential series for the coefficient functions in(x): ∞

in(x) ) exp(n(1 - R)x)

∑ (-1)

m

cn,m exp((m + 1)x)

(7)

m)0

where

c0,m ) xm + 1

(8a)

cn,m ) xm + 1 + n(1 - R)

m



cn-1,k

(8b)

k)0

Notice that each sequence cn,m is a sum of increasing positive values and that the defining sums (8b) can therefore be calculated without numerical difficulty. An equivalent recursion is

c0,m ) xm + 1 for m g 0

This result follows by first observing i0(x) can be expanded as Reinmuth’s series 3b, giving the formula 8a. For n > 0, the result follows by induction, namely, series 7 for in-1(x) is substituted into eq 6b:

in(x) )

d dx



x

exp((1 - R)u) exp(u) + 1

-∞

cn,0/x1 + n(1 - R) ) cn-1,0 for n g 1 cn,m/xm + 1 + n(1 - R) ) cn-1,m +





(-1)mcn-1,m ×

cn,m-1/xm + n(1 - R) for n g 1, m g 1 (9)

m)0

exp((n - 1)(1 - R)u) du exp((m + 1)u)

d ) dx



∫ ∑ x

xπ(x - u)



(-1)k exp(ku)

-∞

k)0

∑ (-1)

m

cn-1,m ×

m)0

exp(n(1 - R)u) du exp((m + 1)u)

xπ(x - u)

The two series are power series in powers of z ) exp(u), and their product can be rewritten according to the convolution product for power series: ∞



∑a z ∑ b k

k)0



mz

k

m+1

)

m)0

m

∑ (∑a

m+1 m-kbk)z

m)0 k)0

The increasing values for cn,m as m increases imply that each series 7 diverges for x g 0; since the values only increase at an algebraic rate, each series will converge for x < 0. EULER’S TRANSFORMATION Although the coefficient functions in(x) exist for all real values of x, as shown by eq 6b, the series expansions of eq 7 only permit the calculation of in(x) for negative values of x. This limitation can be overcome by Euler’s transformation of power series,9 which converts a power series for an analytic function to another series representation for the same function. In the present context, it is convenient to write Euler’s transformation as

f(z) ) ∞





(-1)mamzm+1 S f (z) )

m)0

∑ (-1)

m

( )

(∆ma0)

m)0

z

m+1

1+z

We then have

in(x) )

d dx



where ∆ma0 is the m th repeated difference:

m

∫ ∑ (-1) ∑c x

m

n-1,k

×

-∞

m)0

k)0

exp(n(1 - R)u) du exp((m + 1)u)



m

∑ (-1) ∑c

d

m

)

m)0

n-1,k

∫ dx

x

xπ(x - u)

×

() (

xπ(x - u) ∞



m

∑c

(-1)m

m)0

k)0

n-1,k

xm + 1 + n(1 - R) × exp((m + 1 + n(1 - R))x)



∑ (-1)



exp(n(1 - R)x)

∑ (-1)

m

(∆mcn,0)(1 + exp(-x))-m-1 (11)

m)0

The repeated differences forming the Euler coefficients can be calculated either from the explicit formula 10 or by an algorithm such as

(1) Let a0, a1, ..., aP be given.

m

)

()

In the present case, we apply Euler’s transformation to the series expansions 7 with z ) exp(x), obtaining the new expansions

exp((m + 1 + n(1 - R))u) du

)

)

m m m a a + ... + (-1)m a (10) m m m - 1 m-1 0 0

in(x) )

-∞

k)0

∆ma0 )

cn,m exp((m + 1 + n(1 - R))x)

m)0

(2) For i ) 0, 1, ..., P: set xi ) ai.

where the last two lines represent the induction step yielding formula 8b for cn,m. Here we have used the integral (λ > 0):



x

-∞

exp(λu) du

xπ(x - u)

)

1 exp(λx) xλ

(3) For i ) 1, 2, ..., P; for j ) P, P - 1, ..., i: set xj ) xj - xj-1. (9) Bromwich, T. J. An Introduction to the Theory of Infinite Series, 2nd ed.; Macmillan: London, 1942; Chapter 3.

Analytical Chemistry, Vol. 69, No. 7, April 1, 1997

1467

(4) The resulting sequence is xi ) ∆ia0 for i ) 0, 1, ..., P. In calculating the Euler coefficients, it is essential to realize that the process of forming repeated differences is numerically unstable, so that the calculated value of ∆ma0 results in the true value plus some oscillatory, exponentially increasing error. The error initially results from the truncation error in the fixed precision representation of floating point numbers in computer arithmetic (nominally 16 decimal digits in double precision); this original error then grows due to the numerical instability of the repeated differencing. If the true values of ∆ma0 tend to 0 as m increases, then they are quickly obliterated by the increasing error. Consequently, implementation of Euler’s transformation in traditional computer arithmetic results in Euler coefficients with an unknown, but rapidly growing, error component, and the procedure is of little practical use. Such instability occurs regardless of whether eq 10 or algorithms such as the one above are used. Modern mathematical software provides a form of computer arithmetic which allows practical implementation of Euler’s transformation. Arbitrary precision arithmetic has two important properties: (1) the precision, i.e., length of the floating point representations, for the input values of a calculation can be arbitrarily specified; (2) the precision of the output of arithmetic operations is automatically adjusted to allow for precision loss. When the Euler coefficients are calculated using arbitrary precision arithmetic, the resulting values are formed only to the extent of the precision possible, and the numerical instability is manifested as a steady decrease in precision as higher differences are formed, that is, as more Euler coefficients are calculated. Consequently, implementation of Euler’s transformation in arbitrary precision arithmetic results in Euler coefficients correct to the extent of their accompanying precision, but with that precision gradually decreasing until the calculation is ended by total loss of precision. More Euler coefficients can, of course, be calculated by starting with higher precision representations for the original coefficients. DESCRIPTION OF CALCULATIONS Given R and Λ, we wish to calculate the curve iQR(x) versus x using eq 5. The calculation should be separated into two parts. The first part calculates the coefficients cn,m, using eq 8 (or eq 9), and the corresponding Euler coefficients en,m ) ∆mcn,0, using the algorithm following eq 11. This calculation has the advantage of depending only on R: having done this calculation once, the Euler coefficients can be stored and used with different values of Λ. Arbitrary precision arithmetic is required. The second part calculates the functions in(x) using the Euler coefficients en,m ) ∆mcn,0 in eq 11; the in(x) are then used in eq 5 to obtain iQR(x). This calculation can be done in ordinary double-precision arithmetic. We illustrate the first and second parts with sample calculations on a 33 MHz PC using the general purpose mathematical software MATHEMATICA (Version 2.2.2); in particular, its arbitrary precision arithmetic is used for the first part. (1a) Pick R. The calculations of this paper use R ) 0.3, 0.5, and 0.7. Computational time and behavior were practically the same for all three cases. (1b) Calculate cn,m for 0 e n e N - 1 and 0 e m e M - 1 using arbitrary precision arithmetic. The computation times for 1468

Analytical Chemistry, Vol. 69, No. 7, April 1, 1997

Figure 1. Coefficient functions in(x) for R ) 0.5 calculated by the Euler transformation: largest, i0(x) ) iRV(x); intermediate, i1(x); smallest, i2(x).

N ) 10 and M ) 60, 120, and 180 with 100-digit precision were 15, 30, and 47 s, respectively (R ) 0.5). No loss in precision occurs because the coefficients are generated by adding positive values. (1c) Calculate the Euler coefficients en,m ) ∆mcn,0 for 0 e n e N - 1 and 0 e m e M - 1 using arbitrary precision arithmetic. The computation times for N ) 10 and M ) 60, 120, and 180 with 100-digit precision were 31, 138, and 350 s, respectively (R ) 0.5). As discussed above, substantial loss of precision occurs. The worst case is n ) 9. Using 100-digit precision for the original coefficients in step 1b, the resulting precisions in e9,m were 56, 30, and 8 decimal digits for m + 1 ) 60, 120, and 180. (2) The coefficient functions in(x) were calculated for 0 e n e N - 1 using M-1

in(x) = exp(n(1 - R)x)

∑ (-1)

m

en,m(1 + exp(-x))-m-1

m)0

(12)

For a fixed M, this approximation becomes less accurate as x increases. For example, empirical examination of this approximation at a fixed value of x as the number of terms in the sum changes indicates that, at M ) 180, the absolute errors in i0(x), i3(x), i6(x), and i9(x) are 10-6, 10-9, 10-10, and 10-11, respectively, at x ) 3 and 10-4, 10-7, 10-8, and 10-9, respectively, at x ) 4. Arbitrary precision was used for these accuracy checks, but ordinary double precision can be used in practice. For the approximation N-1

iQR(x) )

∑ (-1) Λ n

-n

in(x)

(13)

n)0

with N ) 10, M ) 180, the difference between arbitrary precision and double-precision calculations over -∞ < x e 4 is extremely small for negative x and increases to about 10-4 near x ) 3 and about 3 × 10-4 near x ) 4. Using double precision and N ) 10, M ) 180, the time to generate one iQR(x) value is about 3.5 s. Taking M ) 180 terms is necessary near x ) 4; the number of terms, and consequently the computation time, can be substantially reduced for smaller x. RESULTS Figure 1 shows the initial coefficient functions i0(x), i1(x), and

Table 1. Sample Values of iQR(x) for Λ ) 2 Using Eqs 12 and 13 with N ) 10, M ) 180 iQR(x)

Table 2. Comparison of Peak Coordinates for r ) 0.5 Λ 2

x

R ) 0.3

R ) 0.5

R ) 0.7

-2 -1 0 1 2 3 4

0.1003 0.1991 0.3156 0.3902 0.3964 0.3621 0.3200

0.0961 0.1952 0.3217 0.4099 0.4173 0.3747 0.3248

0.0908 0.1912 0.3286 0.4275 0.4306 0.3789 0.3240

5

E1/2 - Epeak (mV) ipeak E1/2 - Epeak (mV) ipeak

a

b

c

40.5 0.4232 33.4 0.4361

40.6 0.4232 33.4 0.4361

40.5 0.423 33.4 0.436

a Equations 12 and 13 with N ) 10, M ) 180. b Myland and Oldham series (from Table 2, ref 7). c DigiSim Software (from Table 2, ref 7).

i2(x). Since the coefficient i0(x) is the same as the reversible current iRV(x), it is positive for all x. However, the higher coefficients in(x) show both positive and negative values. Matsuda and Ayabe4 conjectured, but without specific discussion, that the coefficients satisfied a bound

|in(x)/i0(x)| e 0.43 which would imply convergence of expansion 5 for 1 < Λ < ∞. The present calculations support such a bound, and we find for n ) 1, 2, ..., 9 the respective bounds

|in(x)/i0(x)| e 0.431, 0.192, 0.089, 0.041, 0.020, 0.020, 0.011, 0.011, 0.011 Although the calculations support the convergence of the expansion 5, it must be emphasized that this result is an empirical observation based on the first few terms. The possibility that eq 5 is an asymptotic expansion with respect to Λ f ∞ should be considered. In this case, Pade´ approximants7 or a nonlinear sequence transformation10 may be useful in evaluating the series. By empirical examination of the convergence behavior, we estimate that N ) 10, M ) 180 in eqs 12 and 13 yields iQR(x) over the ranges 2 e Λ < ∞, 0.3 e R e 0.7, and -∞ < x e 4 with absolute error on the order of 10-4. Table 1 provides a sample selection of values. Table 2 compares sample values from this paper with those of ref 7, which contains results from both the (10) Basha, C. A.; Sangaranarayanan, M. V. J. Electroanal. Chem. 1989, 261, 431-436. (11) (a) Rudolph, M.; Reddy, D. P.; Feldberg, S. W. DigiSim 1.0 Software; BioAnalytical Systems Inc.: West Lafayette, IN, 1993. (b) Rudolph, M.; Reddy, D. P.; Feldberg, S. W. Anal. Chem. 1994, 66, 589A-600A.

Figure 2. Curve on left: reversible current iRV(x). Cluster of curves: iQR(x) for Λ ) 2; R ) 0.3 (top), R ) 0.5 (middle), and R ) 0.7 (bottom).

Myland and Oldham series and the DigiSim Software11 for electrochemical calculations. Although the work of ref 7 was formulated for a reduction process, the values can be directly translated for comparison here. The conversion from x to E1/2 E uses the standard value f ) 38.9215 V-1 and ne ) 1 in eq 1. Figure 2 shows examples of the results: three curves for quasireversible linear potential sweep current iQR(x) for R ) 0.3, 0.5, and 0.7 (Λ ) 2). The reversible linear potential sweep current iRV(x) is shown for comparison (Λ ) ∞). ACKNOWLEDGMENT The author gratefully acknowledges the helpful comments of Professor Keith Oldham. Received for review August 14, 1996. Accepted December 30, 1996.X AC960826C X

Abstract published in Advance ACS Abstracts, February 15, 1997.

Analytical Chemistry, Vol. 69, No. 7, April 1, 1997

1469