Sensitivity analysis of oscillating reactions. 1. The period of the

A method is developed for the computation of the dependence of the period of an ... may relate the sensitivity coefficients of the period to those of ...
0 downloads 0 Views 521KB Size
J. Phys. Chem. 1981, 85, 1555-1558

In ethylamine and n-propylamine for which the isolated pair spur model is valid, the t values of the initial Gaussian distance distribution have been obtained by both fittings of N ( t ) and F ( S ) and are 40 f 10 A and 50 f 10 A, respectively. From a Dirac initial distance distribution function centered at ro and from the experimental value of N ( m ) , 2 the Rice and Baird relation14yields 27 and 40 A, respectively. Both sets of values may be compared, keeping in mind that the discrepancy is due to the methods of their calculation. In these amines of low polarity, 7 can be also compared to the thermalization range of the solvated electron found to be about 40 A in ether^.^^,^^ By contrast, in water t is either 2340or 27 A41depending upon (38) J. P. Dodelet and G. R. Freeman, Can. J. Chem., 53, 1263 (1975). (39) F. Billiau, J. Belloni, J. A. Delaire, and M. 0. Delcourt, J.Chim. Phys., 76, 1059 (1979). (40) H. A. Schwarz, J. Chem. Phys., 55,3647 (1971). (41) A. Mozumder, J. Chem. Phys., 50, 3153 (1969). (42) W. A. Seddon, J. W. Fletcher, and F. C. Sopchyshyn, Can. J. Chem., 56, 839 (1978).

1555

the model chosen. Thus the r found for ethylamine and n-propylamine conform with the general idea that, the lower the polarity of the medium, the higher the mean distance between initial ion pairs. To the extent that the experimental conditions give isolated pairs, the numerical solution of the one-pair Smoluchowski equation fits well with the experimental data as shown for ethylamine and n-propylamine. For liquids for which the electron undergoes a competitive reaction with the radical, a new model has to be chosen. This work is in progress.

Acknowledgment. This work has been partially supported by the financial specific program A.T.P. 2988 under the auspices of CNRS. The authors are pleased to thank Dr. P. Penel and Dr. A. M. Pelissier of the State University of Toulon for helpful discussions of the mathematical part. (43) W. A. Seddon, private communication. (44) W. A. Seddon, J. W. Fletcher, and F. C. Sopchyshyn, Can. J. Chem.,54, 2807 (1976).

Sensitivity Analysis of Oscillating Reactions. 1. The Period of the Oregonator David Edelson” and Valerie M. Thomast Bell Laboratories, Murray Hill, New Jersey 07974 (Received: December 8, 1980; In Final Form: February 20, 1981)

A method is developed for the computation of the dependence of the period of an oscillating reaction mechanism upon the input parameters of the kinetic model. By expressing the solution of the mass-action differential equations as a delay equation, one may relate the sensitivity coefficients of the period to those of the chemical species concentrations; the latter can be computed by any of a number of established techniques. Application is made to the five-step Oregonator model proposed by Field and Noyes to represent the essential characteristics of the complex Belousov-Zhabotinsky reaction. The method is validated by comparing the results with those obtainable by direct parameter variation. The calculated period sensitivities confirm quantitatively the relative significance of the various components of the model which previous workers could only suggest on the basis of approximate analyses.

The mechanisms which have thus far been conceived to account for the behavior of real oscillating chemical systems,l such as the Bray-Liebhafsky,2 Belousov-Zhabotin~ky,~ and - ~ other bromate systems6J are composed of reaction steps for which the level of independent supporting information has the widest possible range. Some reactions, primarily the inorganic ones, have been measured independently so that their rate parameters are assigned with a fair degree of certainty. Others are estimated from thermodynamic data. The basis of the organic reactions, however, ranges from reasonable confidence to pure speculation, and even the participation of some of the proposed radical intermediates is not assured. Although the mechanisms manage to simulate the gross features of the reacting systems, significant discrepancies exist in detail; furthermore, the uniqueness of the assignment always remains in doubt. The interactions themselves are so complex that the kineticist’s traditional intuition is unsuccessful in placing responsibility upon particular components of the mechanism for an observed feature of the experiment. f

Swarthmore College, Swarthmore, PA 19081. 0022-3654/81/2085-1555$01.25/0

In the past decade techniques for the study of complex chemical reactions by mathematical modeling have advanced rapidly.* These methods are mainly based on the solution of the mass-action differential equations derivable from the chemical model and also provide a route for the study of the effect of the parameters of the model upon its behavior. The use of this “sensitivity analysis” was proposed many years agogand has been well developed for equations capable of analytic solution.1° Numerical extensions with particular application to chemistry have recently been advanced.l1-l3 These could be used to gain (1) Noyes, R. M. Ber. Bunsenges. Phys. Chem. 1980,84, 295. (2) Edelson, D.; Noyes, R. M. J. Phys. Chem. 1979, 83, 212. (3) Field, R. J.; Koros, E.; Noyes, R. M. J.Am. Chem. SOC.1972,94, 8649. (4) Edelson, D.; Field, R. J.; Noyes, R. M. Int. J. Chem. Kinet. 1975, 7, 417. (5) Edelson, D.; Noyes, R. M.; Field, R. J. Znt. J. Chem. Kinet. 1979, 11, 155. (6) Noyes, R. M. J.Am. Chem. SOC. 1980,102, 4644. (7) Orbin, M.; Koros, E. J.Phys. Chem. 1978,82, 1672. (8) Edelson, D. J. Chem. Educ. 1975,52, 642. (9) Poincar6, H. “New Methods of Celestial Mechanics”, 1892-1899; NASA Technical Translation TT F-450-2, 1967. (10) Hille, E. ”Lectures on Ordinary Differential Equations”; Addison-Wesley: Menlo Park,CA, 1969; Chapter 3.

0 1981 American Chemical Society

1556

The Journal of Physical Chemistty, Vol. 85,

Edelson and Thomas

No. 11, 1981

TABLE I: Oregonator Model Parameters’ R a t e Constants h, = 4.00 x 107

h , = 1.34 h , = 1.60 x 10’ h , = 8.00 x 103

where [nIi is a component of the solution exhibiting oscillatory behavior of period T , and the generalized set of parameters a includes the initial conditions as well as the rate constants. If eq 2 is differentiated with respect to the kth parameter CUk

h , = 1.00

Initial Conditions

[ A ] ,= 6.00 [ X I , = 5.025 x a

lo-”

[YI, = 3.00 x 10-7

d[n]i(t,a)/dak = d[n]i(t+~,a) /dak

[ Z ] , = 2.41 X

(3)

T i m e a n d concentrations i n arbitrary units.

further insight into the complex interactions taking place as well as for delineating those components of the model for which additional independent supporting information would be required. The cost of performing these computations for even moderately sized chemical mechanisms has been so large, however, that the technique had been little used. Some newer mathematical improvements’* and the advent of very high speed vector machines (“~upercomputers”)~~ have changed this situation, and notable successes have been achieved in recent applications.16 The usual mathematical formulation of sensitivity analysis for a chemical mechanism seeks the dependence of the concentration of a particular reactant or product at a specified time upon the input parameters of the system, i.e., the rate constants and the initial conditions. In the case of oscillatory reactions, however, it is the characteristics of the oscillations which are of interest rather than specific values of the species concentrations. The stability of the system, i.e., whether or not it does oscillate, the period of the oscillation, the nature of the induction period if any, and the shape of the limit cycle are typical features of interest. It is not apparent how the dependence of any of these factors upon system parameters may be established, since no measure of them appears explicitly in the differential equations describing the system. In this paper we devise a method for relating the sensitivity of one of these factors, the period, to species sensitivities which are obtainable by the usual analysis. The “Oregonator”, a small oscillating reaction system devised by Noyes and Field” to mimic the Belousov-Zhabotinsky reaction, is used to test the method. The parameters which determine the period of this system are established, and the foundation is laid for future application to more complex oscillating mechanisms.

Method of Calculation Assume a mechanism for a chemical reaction system of r reactions in n species formulated in the usual mass-action manner as a set of ordinary differential equations (ODES): r

d[n]i/dt = C v . a j=’

rJ

J

f [n]lyI)

1=1

i = 1, 2, ..., n

(1)

where [n]l is the concentration of the Zth species, aj is the rate constant for the j t h reaction, and vij is the stoichiometric coefficient of the ith species in the jth reaction. If the mechanism is that of a chemical oscillator, then the solution of this set of ODES satisfies the delay equation (11)Atherton, R. W.;Schainker, R. B.; Ducot, E. R. AIChE J. 1976,

21, 441.

(12)Cukier, R. I.; Levine, H. B.; Shuler, K. E. J.Comput. Phys. 1978, 26,1 and references cited therein. (13)Dickinson, R. P.;Gelinas, R. J. J. Comput. Phys. 1976,21, 123. (14)Hwang, J.-T.;Dougherty, E. P.; Rabitz, S.; Rabitz, H. J. Chem. Phys. 1978,69, 5180. (15) Edelson, D. ACS Symp. Ser., in press. (16)Edelson, D.; Allara, D. L. Int. J. Chem. Kinet. 1980, 12, 605. (17)Field, R. J.; Noyes, R. M. J . Chem. Phys. 1974, 60, 1877.

relating the desired sensitivity of the period T with respect to parameter a k to the usual species sensitivity a[nIi/aak. The computational procedure to be used now becomes apparent: (1)The mass-action ODES are solved for conditions which exhibit reproducible oscillations. (2) A suitable oscillating species is selected and its sensitivity to the parameters f f k evaluated at corresponding time points in two adjacent cycles. (3) The period sensitivity is evaluated from eq 4.

Oscillating System The Oregonator of Noyes and Field was used for the purpose of the present study. This five-reaction system was devised to exhibit the behavior of the three main processes of the Belousov-Zhabotinsky system without getting involved in the detailed chemistry. The model, in one of its variations, is A+Y-X

(R1)

X+Y-P A X -2X

032)

+

2X

-

z-Y

Q

+Z

(R3) 034) (R5)

With the rate constants and the initial conditions of Table I, the three intermediates X, Y, and Z exhibit oscillatory behavior of Figure 1. The starting material A is gradually depleted, while the end products P and Q accumulate. A closed system will exhibit damped oscillations with a gradually lengthening period. The simulation was therefore performed for an open system in which A is kept nearly constant by replenishment to its initial value. This was accomplished by adding a term of the form k,( [A], [A]) to the differential equation for [A]. In addition, small amounts of X, Y, and Z were introduced at the beginning to eliminate any induction period, as indicated in Table I. As a result the oscillations showed a highly reproducible wave shape with a period of ca. 9 time units. In principle any of the oscillating intermediates would be suitable for performing the sensitivity analysis. However, when the form of eq 4 is considered, it is apparent that X would be a poor candidate, since the low slope of the concentration curve between the autocatalytic spikes would lead to poor accuracy for the difference in the numerator, and this error would be further magnified by division by the small quantity in the denominator. Y would be better, but the large negative spikes could be

The Journal of Physical Chemistty, Vol. .85, No. 11, 198 1

Sensitivity Analysis of Oscillating Reactions

1557

TABLE 11: Oregonator Period Sensitivity Coefficients a= aT/aa

A T / ( A ~ 0.5%) = AT/(A= ~ 1.0%) A T / ( A =~ 2.0%)

-8.104 -8.527 -8.324 -8.338

x

x lo-' 5.823 x x lo-' 6.263 x x 10-' 6.077 X x lo-' 6.128 x

lo-'' lo-'' lo-'' lo-''

Period 2 3.130 x 6.106 X 5.705 X 5.687 X

x x

a In r/a In a

-1.204 X

5.967 6.207 6.337 6.166

lo-' 10-l lo-'

lo-'

X X X

1.058 X 10.'

'

4.155 X

I

1

k,

lo-'' lo-'' lo-'' lo-''

x lo-' X

A T / ( A=~ 0.5%) A T / ( A ~ 1.0%) = A ~ / ( A = ~2.0%)

-9

k,

Period 1 -1.663 X 4.683 X 6.325 X l o M 5 -1.866 X -1.810 X 6.034 X 5.808 X l o M 5 -1.838 X

-7,908 -8.670 -8.320 -8.442

aT/aa

k'

kl

I

I

I

I

1

TIME

Flgure 1. Oscillations of the concentrations of the Oregonator intermediates.

troublesome. Z, on the other hand, appears to be the best behaved and was chosen as the variable to be used in the analysis.

Computation Method The species sensitivities required in eq 4 were computed by using the Greens function method (GFM) of Rabitz et aL1* implemented as previously described15on a Cray-1 computer (although any of the other established sensitivity methods could be used). Because of the extremely steep edge and sharp corners of the concentration wave form, it was found necessary to perform the ODE solution for the Greens function with very stringent accuracy requirements. Furthermore, the trapezoidal rule integration used does not perform well at the almost delta function shape of the integrand at the oscillator switch points. While it should be possible to obtain comparable results by using species sensitivities at any two adjacent cycles, this factor led to a progressively more serious degradation

-1.480 x

lo-* lo-' lo-' lo-'

-8.921 -8.904 -8.882

lo-'

-6.498 X -6.483 X -6.558 X

-8.810

-1.825 X

lo-'

-8.381 -8.948 -8.891 -8.814

-7.377 X

lo-'

-9.894 X

-1.808 X l o - ' -1.811 X 10.'

lo-'

[AI,

k,

-6.356 X -6.767 X -6.733 X

lo-'

-4.324 X

lo-' lo-' lo-' lo-' lo-'

lo-'

of accuracy as the integration proceeded through successive cycles, and it was found necessary to limit the computation to the first few cycles. In addition to the analysis by the GFM, sensitivities of the period to input parameters were obtained by direct variation of the rate constants and initial conditions, as well as by an intermediate scheme in which species sensitivities obtained by direct variation were used in eq 4 to give the period sensitivity. These served as comparison calculations to validate the formal analysis. Since the GFM yields the first-order sensitivity coefficients, the direct variation was made for a number of values of Acu so that the effect of higher order terms could be eliminated in making the comparison.

Results Table I1 lists the sensitivity values dT/dcu obtained by the above analysis, for the parameter values in Table I, where CY represents any of the five rate constants or the initial condition [A],. The evaluations were made over the first (t = 4.367,13.383) and the second ( t = 13.383, 22.400) periods. For comparison, sensitivities obtained by direct variation of cu (by the indicated percentage of original value) are also given. As expected, the direct sensitivities are in good agreement from one period to another, and for the most part either scatter about or tend toward the computed value for period 1. The computed value for period 2 is in poorer agreement, and this degradation in accuracy gets worse as the analysis is extended to longer times, owing to the progressive accumulation of numerical error in integrating the Greens function through successive cycles. As a further test of the technique, sensitivities calculated by eq 4 from A[Z]/Acu values obtained by direct variation were in good agreement with those in Table 11. The final line in Table I1 expresses the period 1computed sensitivitiesin logarithmic or relative (dimensionless)form. Discussion The results obtained in this work indicate that the principle factor determining the period of oscillation of the Oregonator is k6,the rate of reaction R5 which corresponds to the bromide regeneration process of the BelousovZhabotinsky mechanism. This reaction is almost an order of magnitude more significant than the others, with reactions R3 and R4 being of least importance. The initial concentration of A is also of secondary importance. Although the Oregonator model has been extensively analyzed mathematically,1s no one has derived an expression for the period in terms of the systems parameters. Field and Noyeslgqualitatively estimate that the "delayed feedback" represented by reaction R5 will be the main determining factor and that reactions R3 and R4 will be (18)Hastings, S. P.; Murray, J. D. SIAM J . Appl. Math. 1976,28,678. (19) Field, R. J.; Noyes, R. M. Acc. Chem. Res. 1977, 10, 214.

J. Phys. Chem. 1981, 85, 1558-1562

1558

the least important since they are components of a fast switch which is but a small part of the total period. In a more elaborate analysis, Tyson20finds that the length of parts of the cycle are related to the stoichiometric factor h (Z hY; here h = 1)in reaction R5 again emphasizing this reaction as the principal determinant of the oscillation period. These earlier qualitative considerations may now be considered to be confirmed by the results of the present paper.

-

(20) Tyson, J. J. “Lecture Notes in Biomathematics”;Springer-Verlag: New York, 1976; Vol. 10.

The other main result of this work is the validation of an easily implemented technique for studying the factors determining the period in the much more complicated real systems of which the Oregonator is a skeletal model. In the Belousov-Zhabotinsky system, the regeneration mechanism involves several organic reactions which are highly speculative, and it is anticipated that the application of sensitivity analysis to that mechanism will greatly assist in the clarification of its remaining unresolved features. Acknowledgment. We are indebted to N. L. Schryer for suggesting the delay equation formulation of the sensitivity analysis used in this paper.

Photodissociation Spectroscopy of Gas-Phase Protonated Acrolein Derivatives Jeffrey P. Honovlch and Robert C. Dunbar”? Department of Chemistry, Case Western Reserve University, Cleveland, Ohio 44 106 (Received: December 29, 1980)

The low-resolutionphotodissociation spectra of a series of a,b-unsaturated aldehydes and ketones have been obtained by using ion cyclotron resonance techniques. Comparison of these spectra with the absorption spectra of the corresponding neutrals yields the acid-base properties of the aa* excited state in terms of a change in proton affinity (APA) upon excitation. The results concern the effects of methyl groups and extended conjugation on the proton affinity of the *aa*state of the parent acrolein chromophore. Some of the ions in this study are compared with solution spectra to determine the effects of solvation. The results of several theoretical calculations are examined in terms of excited-state properties and compared with our APA interpretations. The noted similarities and differences between ion and neutral indicate that a set of modified Woodward-Fieser rules would be appropriate for this class of ions.

Introduction Compared to the physical and chemical knowledge of neutral species, our understanding of the corresponding ions is far less complete. Advances in experiment and theory are now allowing many of the approaches worked out for characterizing neutrals to be applied to these less familiar compounds. Of particular interest to us has been the recent emergence of optical spectroscopic means of studying gas-phase ions made possible by the technique of ion cyclotron resonance photodissociation spectroscopy (ICR-PDS). We describe here a study of a series of protonated a,P-unsaturated aldehydes and ketones by this technique where the results have connections to the well-established spectroscopy of the corresponding neutrals. The influence of substituents on the ?r-a* transitions in neutral chromophores is documented in large detail, with some of the regular and systematic effects being codified in the empirical Woodward ru1es.l As yet, no sufficiently extensive series of gas-phase ions has been examined to assess how these quantitative and qualitative regularities carry over to the ions. The protonated acrolein chromophoric group was chosen for such a study, with the hope that protonation of the oxygen would not disturb the neutral chromophore beyond recognition, so that comparison with the well-understood neutrals would be possible as a reference. Experimental Section Due to the difficulty in obtaining large enough concentrations of gas-phase ions for normal absorption spect John

Simon Guggenheim Memorial Fellow, 1978-79. 0022-365418112085-1558$01.25/0

troscopy, the indirect method of photodissociation spectroscopy (PDS) was employed as the route to obtaining spectroscopic information. As details of ICR spectroscopy are given elsewhere,2 only a brief description will be presented below. The ICR instrument represents a specialized mass spectrometer which functions to generate ions, trap them on a time scale of seconds via crossed electric and magnetic fields, and finally detect them according to their m l z ratio. During the trapping period the ions are irradiated with wavelength-selected light. If absorption occurs and the ion reaches either directly or indirectly a dissociative state, the result will be a decrease in the population of the absorbing ions. By monitoring the disappearance of a specific ion and normalizing for irradiation time and intensity, one obtains a photodissociation cross section (a), which, when plotted as a function of irradiation wavelength, yields the photodissociation spectrum of the ion. Anything observed by this method represents at most the lower limit of the absorption spectrum of the ion, since photodissociation may be in competition with nondissociative relaxation processes such as fluorescence. The quantum yield for photodissociation is usually taken to be one or a constant over the absorption band. As long as the photon energy is well above the thermochemical dissociation threshold (as is the case here, where loss of H20 or H2 can be estimated to be less than 2-eV endothermic),there is no reason to anticipate a strong wavelength dependence of the dissociation quantum yiel.d, nor have the large variety of polyatomic ions studied by (1)A. I. Scott, “Interpretation of the Ultraviolet Spectra of Natural Products”, Pergamon Press, Elmsford, N.Y., 1964. (2) R. C. Dunbar, E. W. Fu, and G. A. Olah,J.Am. Chem. SOC.,99, 7502 (1977).

0 1981 American Chemical Society