Computer construction of reaction mechanisms - The Journal of

Computer construction of reaction mechanisms. Benny G. Johnson, and P. L. Corio. J. Phys. Chem. , 1993, 97 (46), pp 12100–12105. DOI: 10.1021/ ...
1 downloads 0 Views 669KB Size
J. Phys. Chem. 1993,97, 12100-12105

12100

Computer Construction of Reaction Mechanisms Benny G. Johnson+and P. L. Corio' Department of Chemistry, University of Kentucky, Lexington, Kentucky 40506-0055 Received: June 24, 1993; In Final Form: September 3, 1993'

A method is given for the computer construction of all possible reaction mechanisms without reentrant steps. The maximum molecularity may be set a t any desired value, and all possible mechanisms consistent with this choice are then constructed. The reaction components must be specified, but the procedure is otherwise independent of the energy and physical states of the components. It is assumed that the reaction occurs in a closed system, but as no assumption is made regarding the energy, the reaction need not occur isothermally.

Introduction Reaction mechanisms have been a cornerstone of chemistry for almost a century and, for the most part, have been constructed by considerations based on energetics, kinetics, and intuition. Although such considerations are important, there is always the question as to the existence of rival mechanisms that may have been overlooked. Nevertheless, the problem of constructing all possible reaction mechanisms consistent with a given set of components has only recently been addressed,Iv2even though there has been some progress regarding the structure of reaction mechanisms.)" It has been shown192 that reaction mechanisms are subject to certain equalities and inequalities that are useful in the construction problem. In particular, K = P + T + L = ~ + v

(1)

where K is the number of reaction components, p , r, and L, respectively, denote the number of reactants, products, and intermediates, and n and v, respectively, denote the number of independent reactions and number of independent conservation conditions.7 For example, the mechanism proposed for the gasphase decomposition of ozone, namely, O)+O,+O 0,

+0

--*

20,

satisfies these relations, since K = 3, p = r = 1 = v = 1, and n = 2. The values of p and .rr can also be obtained from the stoichiometric relation

20,

-

30,

(3)

but in the case of a catalyzed reaction it is necessary to include all catalysts238 when computing p r. For brevity, mechanisms conforming to eqs 1 and their implications will be called consistent mechanisms. We discuss here theapplication of this theory to thedevelopment of algorithms for generating all possible consistent mechanisms for a given reaction. We have implemented these algorithms in computer programs, written in TURBO C for the IBM PC family, which fully automate the mechanism generation process. To facilitate the mathematical treatment, we shall denote a general chemical reaction

+

(4) j= 1

where the Xj represent the K distinct components and the uj their integral coefficients, after transposing all reactants to the right member. Thus, uj will be positive for a given reaction if Xi is a Present address: Department of Chemistry, Carnegie Mellon University, Pittsburgh, PA 15213. * Abstract published in Advance ACS Absrracrs, October 15, 1993.

0022-3654/93/2097-12100$04.00/0

product, negative if it is a reactant, and zero if its (net) coefficient is zero. When the Xj have been given, the reaction is completely specified by the coefficients u,, and as a concise notation for (4) we introduce the column vector x = (ul,..., where T denotes matrix transposition. The conservation conditions for the system may be written as

where the square brackets denote molar concentration, and the C, are constants determined by the initial conditions. One conservation condition can be written for each chemical element in the system, plus one for conservation of charge, if necessary. This may yield more than v conservation conditions (in which case these form a linearly dependent set of rank v), but we need concern ourselves only with a set of v independent conditions. The ai, are numerical coefficients, and it is always possible to write the system ( 5 ) with all ai, integers; we shall henceforth assume the afj are integral. Let A = (aij)vx,be the matrix of coefficients in ( 5 ) , and let x be a vector representing a mechanistic step of a consistent mechanism. Then it can be shown2 that

Ax = 0 (6) Le., x is in the null space of A. Since the rank of A is v, the dimensionality of its null space is n, and therefore if X = (u,+,, is a matrix of coefficients correspondingto a consistent mechanism, then the columns of X (the individual steps) comprise a basis for the null space of A. The problem of constructing all consistent mechanisms may now be separated into two parts. First, the entire set of valid mechanistic steps is generated, using (6) as a starting point, without considering how to assemble these steps into mechanisms. After this list of steps is obtained, all consistent mechanisms are subsets of it, and at this point the question as to which subsets are appropriate is addressed. Generation of Mechanistic Steps Equation 6 formally admits infinitely many solutions; however, not every solution is physically meaningful. For example, since the components of x are the coefficients in a mechanistic step, only solutions for which these coefficients are integers are admissible. Furthermore, it is reasonable to assume that no step will exceed a certain maximum molecularity, say, two or three. When the maximum molecularity is restricted to some finite value, the solutions of (6) and, hence all possible mechanisms, are finite in number. As every solution of ( 6 ) satisfying these criteria represents a possible mechanistic step, it is necessary to obtain all such solutions in order to ensure that all relevant mechanistic possibilities may be determined. Any such solution may be expressed as a linear 0 1993 American Chemical Society

The Journal of Physical Chemistry, Vol. 97, No. 46, 1993 12101

Computer Construction of Reaction Mechanisms combination of the vectors in any basis for the null space of the matrix A, so that, in principle, all solutions are determined simply by constructing such a basis. For this reason, when generating the desired solutions it will be convenient to work with their coordinates with respect to some basis rather than with the actual components of the solutions themselves. However, it is not feasible, in general, to construct all solutions of the proper form from an arbitrary basis, for even though the components of the basis vectors may consist entirely of integral components, their coordinates are not necessarily integral, or even rational. To make the problem tractable, it essential to obtain a basis for which all solutions of (6) in integers can be constructed from integral linear combinations of the basis vectors. The existence of such a basis is provided by the following theorems of H. J. S . Smith.9 Theorem 1. The System (6) has SOlUtiOnS X k = ((Tlk, ..., u,k)T, k = 1, ..., n, such that the determinants formed from the n X n submatrices of the fundamental matrix X = ( U j k ) have greatest common divisor unity. Theorem 2. These solutions x k form a system of fundamental solutions in the sense that every solution in integers is expressible in the form n

(7) where the Ck are integers. Skolem’s proof10 of these theorems suggests a method for constructing a fundamental matrix solution with the desired properties. Such a matrix will henceforth be called a Smith matrix and the basis comprised by its columns a Smith basis.” Smith’s theorems guarantee that admissible mechanistic steps canbemade to haveintegral coordinates with respect toa specially chosen basis but give no indication of what the values of these coordinates will be. However, the additional criteria that the solutions of (6) satisfy impose certain constraints upon thevalues of the components of the solution vectors. In particular, the maximum molecularity constraint limits the possible forms of the reactant member of a mechanistic step. Let m be the maximum allowed molecularity. Then, for a system with K components, the number r of distinct sets of reactants is given by

r = -( K-

+ m)! K!m!

1

For each of these r reactant sets there may, in general, be zero, one, or more than one combination of the remaining components which, when placed on the products side, yields a valid mechanistic step. The algorithm used for generating all appropriate solutions does so by considering each distinct set of reactants separately and then determining all solutions corresponding to that set. Note particularly that here the term “reactants” is applied only to a particular set of mechanistic steps, namely, those that appear on the left side of certain (elementary) reactions, without regard to their function in the overall reaction, that is, the stoichiometry. Given a particular choice of reactants, all components of the vector x representing any step having those reactants are bounded. The component of each reactant in x must have a definite negative value from -1 to -ma All other species must appear on the right side of that step or not at all, hence their vector components must be greater than or equal to zero. Let R be the current set of reactants, and x = (ai,...,u ~be)a solution ~ of ( 6 ) with the given reactants (note that in some cases no such x exists). Defining

on the vector components. The same number of atoms of each distinct element in the system must appear in each member of a reaction; in particular, the weaker observation that no more mass may appear on the right than on the left can be used to obtain a bound. Let e be the number of distinct elements present in the system and let E,, i = 1, ..., c be those elements. Define eij, i = 1, ..., e, j = 1, ..., K as the number of E, in chemical component Xp Then an upper bound to a j is

,XER a

i

1, ..., K (10) where the incomplete (floor) brackets denote the greatest integer function. Taking (9) and (10) together, we have j = 1, ..., K

ai I ai I a;,

(1 1) with strict equality whenXj is in the reactant set. Thus all vector components can either be fixed or bounded with a definite range. Let k = 1, ..., n be the coordinates in a Smith basis X = ( U j k ) of a solution x whosecomponentslie within thegiven bounds. Any solution can be expressed as a linear combination of the columns of X;therefore, on substituting for x its expansion in the Smith basis, bounds are obtained on linear combinations of the n

Note that, in general, the bounds imposed on the individual Ck are not immediately discernible; to determine them, the system of linear inequalities (12) would have to be solved, a nontrivial task. Fortunately, however, this can be avoided; for the Smith fundamental matrix is not unique, and the following corollaries to Smith’s theorems show how this can be used to one’s advantage. Corollary 1. There exists a Smith matrix for the system (6) which contains, perhaps after permutation of rows, an n X n submatrix L = ( k j k ) , where L is lower triangular. Corollary 2. There exists a Smith matrix as described in Corollary 1 for which k j j > 0 and k j k 5 0, lkjkl < l k j j for k < j, j = 1, ..., n. Constructiveproofs of both corollaries are given in the appendix. Once a Smith matrix X of the type of Corollary 2 is obtained, bounds on the individual coordinates can be readily determined by a straightforward substitution process, Let p = @I, ...,pn) be the permutation vector giving the rows to X which make up the submatrix L, i.e.,

Then the following recursive formulas generate lower and upper bounds c;; and c,: respectively, on the coordinates c k of the desired solutions with respect to X.

k = 2,

...,n I

we see that a; Iu j , i.e., u; is a lower bound to u,. The principle of mass conservation places another constraint

(14) k-1

I

12102 The Journal of Physical Chemistry, Vol. 97, No. 46, 1993

where the incomplete (ceiling) brackets for the lower bounds denote the smallest integer not less than the obtained result. After the coordinate bounds have been obtained from (14), the relevant solutions are easily obtained. All vectors with integral coordinates satisfying the given bounds are constructed. By generating the vectors from their coordinates, the need to check for satisfaction of (6) is avoided; furthermore, Smith's theorems guarantee that no pertinent solution exists which has any nonintegral coordinates. By design, the vectors considered will have components within the ui and u/t for a l l j such that the j t h row of X is in L. It is then necessary only to check that the remaining K - n = Y components fall within the bounds before accepting the solution. Note from (14) that the difference c: - c; for k > 1 depends directly upon the corresponding difference in component bounds :u - u;, the coordinate bound differences for all previous k, as well as the magnitudes of the elements of the submatrix L. The property of a Smith matrix of Corollary 2 that the element of largest magnitude in each row of L lies on the diagonal serves to keep the summations in (14) small and hence the coordinate bound differences are kept small. In practice, the diagonal elements of L are usually all 1 or 2, and thus the differences are usually comparable to the differences between the component bounds used to generate the coordinate bounds. Repeated application of this method is used in a K-stage process to obtain all mechanistic steps possible for a given set of reactants. At the j t h stage, a Smith matrix of Corollary 2 having pI = j is constructed and used as the basis for generating solutions. Of the r choices of reactants for a mechanistic step, those in which Xj appears but XI, ...,Xj-Ido not appear are considered in turn, and all solutions of (6) corresponding to steps having exactly those reactants are generated as prescribed. These stages are mutually exclusive, and it is easily seen that all steps are generated by this process. The motivation for the staged approach is efficiency. We see from (14) that whenpl = j for the Smith matrix at thejth stage,

c; = c1 = c;

(15) Thus, the coordinate ci is uniquely determined, and this serves to minimize the ranges in which the remaining coordinates are bounded. At each new stage, it is not necessary to construct the new Smith matrix by Skolem's method;6 the existing Smith matrix may simply be modified by a Gauss elimination process. As an example, we return to the ozone decomposition mechanism considered earlier. When this algorithm is applied to generate all steps having molecularity 2 or less, the following list is obtained:

0,-O,+O 0,--* 3 0 20, 30, 20,

20,

- ++ + -- + + 20,

20

0,

40

+

20,

60

03+02-+50 0, 0 2 0 , 0,

20,

20

0,

20

0

40

20,

0,

(16)

0

0, 0,

In this case, it can be easily verified by inspection or graphically2

Johnson and Corio that these are indeed all the possibilities. Thevalue 2 was chosen as the maximum molecularity to be practical and also to maintin simplicity of the example. However, the method can accommodate any finite molecularity. Note that some of the steps in (16) are merely the reverses of other steps in the list. Here, a step and its reverse should not be considered together as one step, because the two are distinct solutions of (6) (the reverse solution being the negative of the forward solution). The same is true for steps which are positive multiples of other steps (such as the first and fourth steps in the list), for these represent two quite different events on the microscopic level. These are simply particular examples of a set of steps which are linearly dependent. Since n = 2 here, any set of three or more steps must be linearly dependent. Not all steps in (16) occur, for example, in forward-reverse pairs because the molecularity constraint prohibits the reverse step, even though it does satisfy (6).

Generation of Consistent Mechanisms The fundamental relations (1) tell exactly how many independent steps must appear in any consistent mechanism; however, they give no indication of which combinations of steps should be considered as possible mechanisms. In principle, the elucidation of the entire set of mechanistic steps is the most important part of the problem: since each consistent mechanism is a subset of this set, all possible mechanisms are implicitly determined at this point and could be obtained by taking all subsets of n or more steps, subject to certain constraints which must be satisfied by any mechanism,' e.g., it is obvious that a mechanism which does not incorporate all components is invalid. However, this procedure can yield a vast number of consistent mechanisms for all but the simplest reactions, so at this point, additional assumptions about the mechanisms will be introduced. These are not unique and can vary according to factors such as empirical evidence about the type of mechanism. 1. All of the K components must appear in the mechanism. 2. The mechanism has exactly n steps. 3. The first and last occurrence in the mechanism of all reactants must be on the left side of a mechanistic step. 4. Thefirst andlastoccurrenceinthemechanismofall products must be on the right side of a mechanistic step. 5 . The first occurrence of all intermediates must be on the right side of a mechanistic step, while the last occurrence must be on the left side of a mechanistic step. 6. The first occurrence of all catalysts must be on the left side of a mechanistic step, while the last occurrence must be on the right side of a mechanistic step. These assumptions characterize a nonchain mechanism having no dependent steps. They do not encompass all categories of mechanisms invoked to explain chemical reactions and are not intended to do so. They simply illustrate the types of additional constraints which could be used to delimit the consistent mechanisms selected. The above assumptions are the most restrictive for the very first step of the mechanism (only reactants and catalysts on the left side, only intermediates or products on the right side) and for the very last step (only products and catalysts on the right side, only intermediatesor reactantson theleft side). It is therefore preferable to use a 2-fold approach, starting from both the beginning and end of the mechanism and working inward. All steps satisfying the requirements for the first mechanistic step are selected, and as mentioned above, these criteria depend only upon the functions of the components involved. Then, for each valid first step, all steps from the master list are selected which, when appended to the first step, yield a valid two-step beginning partial mechanism (BPM). In general, determining the second step is more complicated than determining the first. There are usually more possibilities for the second step, for not only many

The Journal of Physical Chemistry, Vol. 97, No. 46, 1993 12103

Computer Construction of Reaction Mechanisms reactants or catalysts appear on the left side but also any intermediates or products generated by the first step; furthermore, the exact possibilities depend on which of the potential first steps is being used. In a similar fashion, ending partial mechanisms (EPMs) are constructed as well. The beginning and ending partial mechanisms are extended by appending steps in this way until the number of beginning plus ending steps reaches n. Representatives from each set of partial mechanisms may then be joined to see if the union gives a mechanism consistent with the theory and additional assumptions. This approach has advantages over simply extending BPMs to n steps; fewer potential mechanisms aregenerated, and no partial mechanism need ever be tested for satisfying the requirements of both a beginning and an ending partial mechanism. To aid in the quantification of this approach into an algorithm’ for generating all mechanisms, we introduce the convention of a status uector of a partial mechanism. This is vector with K integer components corresponding to the chemical components of the system. Let u = (u1, ..., u,) be the status vector associated withagivenBPMandv = ( V I , ...,ux)bethestatusvectorassociated with a given EPM. The componentsof u and v indicate the status of the chemical componentsin the partial mechanisms as follows. For BPMs,

( uj =

6. Check the following conditions for the union of each BPM with each EPM: (a) For n o j E (1, ...,K ) , uj = 0 and uj = 0, and for no XI which is a reactant: uj = 0 and uj = 2, or uj = 2 and uj = 0 product: uj = 0 and uj = 1, or uj = 1 and uj = 0 intermediate: uj = 0 and vi = 2, or uj = 1 and uj = 0 catalyst: uj = 0 and uj = 1, or uj = 2 and uj = 0 (b) The union of the BPM and x is a linearly independent set. If both of these are satisfied, add the union to the list of consistent mechanisms. As an example, when given the set of mechanistic steps (16), this algorithm produces the following list of ten mechanisms for the conversion of ozone to oxygen.

0,-,0,+ 0 O3+O-,2O2

-

0,-O,+O

20

0 e X j does not occur

0,

0, 3 0 -+

1 e X j last occurs on the left

j = 1, ...,K

(2 * Xi last occurs on the right

03+O-.20,

-

0, 3 0 +

and for EPMs,

20

0 =+ X j does not occur uj =

t

1 e X j first occurs on the left

j = 1, ..., K

-+ -.+ -- +

20,

20,

0, 0

2 =+ Xi first occurs on the right (18)

Given a BPM status vector u and a mechanistic step x, the validity of appending x to the BPM can easily be determined by checking the following: 1. For all j E (1, ..., K ) , if Xjis a reactant: uj = 0 and uj < 0 product: uj = 0 and uj > 0 intermediate: uj = 0 and uj > 0 catalyst: uj = 0 and uj < 0

2. The union of the BPM and x is a linearly independent set. If both of the above two conditions are satisfied, the union is itself a valid BPM; otherwise, x may not be appended to the BPM. Analogous criteria are applied for adding a step to an EPM. The following algorithm then generates the set of all consistent mechanisms. 1. From the list of all mechanistic steps, find all which are valid first steps of the mechanism. This is equivalent to applying the BPM extension criteria using a null BPM. 2. Initialize a BPM status vector for each beginning step. 3. From the list of all mechanistic steps, find all which are valid last steps of the mechanism. This is equivalent to applying the EPM extension criteria using a null EPM. 4. Initialize an EPM status vector for each ending step. 5. While the total number of BPM and EPM steps is less than n: If the current number of BPMs is less than the current number of EPMs, create the next set of BPMs by extending the current BPMs by one step as described, selecting the additional steps from the entire list of mechanistic steps, and creating an updated status vector for each new BPM; else, extend the current EPMs analogously.

0,

20,

20

20,

20,

20

0,

20

+

20, -.* 0, 4 0 O3+O+2O,

0,+ 40 2 0 0,

20,

-.*

-

+

20,

60

0, 0 -.* 20, 20, 20

-

-+

60

0,

Since the list contains all consistent mechanisms for the reaction (subject to the additional assumptions), only mechanisms in this list need be considered when searching for the true mechanism. The first mechanismgiven,with a reversible first step, is consistent with experimental2 It is appropriate to note at this point that reversible reactions have not been explicitly taken into account for several reasons. (1) The primary objective of the research was the construction of reaction mechanisms, and the method used to construct them is independent of the reversibility or irreversibility of any mechanistic step, though such considerations would influence the kinetics. (2) The same consideration applies to thegeneration of the list of independent reactions from which the reaction mechanisms are compiled, that is, the same list would be obtained for any assumption as the reversibility or irreversibility of the

Johnson and Corio

12104 The Journal of Physical Chemistry, Vol. 97, No. 46, 1993

reactions. (3) The desire to develop a procedure that, as far as possible, was independent of physical assumptions. In practice, of course, the question of reversibility must be addressed. This could be done formally in a computer program, which would introduce 2" contingencies for each consistent mechanism of n (independent) steps. Alternatively, one could rely on experiment to decide the issue. Theoretically, the principle of microscopic reversibility could be invoked as a guide to reversibility, keeping in mind that the principle applies only to systems in thermodynamic equilibrium, so that some caution must be exercised in its application to systems displaced from equilibrium. Nevertheless, excepting cases where the reverse reaction is of such high molecularity as to render them improbable (e.g., the reverse of the first steps of mechanisms 5-10, above) the principle can be useful in deciding the question of reversibility. Concluding Remarks There are four advantages of employing this approach as an aid to mechanism determination. (1) The set of all consistent mechanisms can be generated without regard to experimental detail (apart from that which may be necessary to determine the chemical species in the system); the method is independent of the energetics involved. Given the reactants, products, catalysts, and known or proposed intermediates of a chemical reaction, a meaningful mechanistic study can be carried out by applying the simple but powerful theory which follows solely from the laws of mass and charge conservation. (2) As can be demonstrated by example,' mechanisms have been presented in the literature which do not satisfy the fundamental relations of the theory, and are therefore inconsistent. Such mechanisms certainly would not have been proposed did they not seem to explain features of the kinetics; this shows that in this way it is possible for an inconsistent mechanism to pass as a valid one. However, if only consistent mechanisms are considered, the possibility of being misled by an ostensibly feasible inconsistent mechanism is automatically excluded. (3) As the true mechanism must itself be consistent, it is guaranteed to be contained in the list of all consistent mechanisms generated by this method (provided that any additional assumptions imposed on the mechanisms are valid). This was evidenced in the ozone example. (4) While consistency does not guarantee that a mechanism will explain the significant features of the kinetics, the list of all consistent mechanisms can be used as a source of potential viable alternatives to a proposed mechanism that might have otherwise been overlooked. A major drawback of the method as currently implemented is that certain of the additional assumptions about the mechanism are nonrigorous. Perhaps the most restrictive of the assumptions listed earlier is that there must be exactly n steps in the mechanism. Therefore, the current method is not able to generate mechanisms having dependent steps, except trivially by allowing the reversibility of individual steps. A notable class of mechanisms excluded by this restriction is that of cyclic mechanisms. Of course, any dependent steps are expressible as linear combinations of the n chosen independent steps. We stress that the limitation on the types of mechanisms which can currently be analyzed in no way results from the manner in which the individual mechanistic steps are determined but is solely a consequence of what follows after the list of all steps is obtained. The step generation algorithm was developed strictly from the conservation laws with no approximation other than finite maximum molecularity, which can always be increased to a value great enough to encompass all physically practical possibilities, and hence is entirely rigorous. The systematic elucidation of all mechanistic steps involving a given set of chemical components is therefore a solved problem. The ozone reaction, while sufficient to illustrate the method, is an extremely simple example. For more complicated reactions, the number of consistent mechanisms will generally be substrin-

tially larger. For example, it would not be impossible for a reaction having only 10 components to have hundreds of possible mechanistic steps and over a million consistent mechanistic possibilities, even when imposing the additional assumptions. This is therefore the proper point to introduce the energetics of the reaction, as a means of delimiting the mechanistic steps and mechanisms obtained from the theory, rather than in an ad hoc constructive approach. The logical extension of this work is briefly outlined as follows. Mechanistic steps which are thermodynamically unfavorable can easily be eliminated by computing approximate reaction enthalpies from a database of experimental heats of formation. Likewise, comparison against a library of experimental activation energies allows steps which are kinetically prohibited to be identified. In the cases where reliable experimental data are not available, quantum chemical ab initio methods can be employed instead to calculate the desired results. This will reduce the total number of mechanistic steps considerably, and will allow nonrigorous assumptions precluding certain mechanistic types to be relaxed without the total number of mechanistic possibilities becoming unmanageable. This will also facilitate the development of more sophisticated algorithms capable of considering, e.g., chain mechanisms and cyclic mechanisms. Once all possible consistent mechanisms are obtained from the pared list of steps, each can be checked against experimental kinetic profiles, using experimentally known or theoretically calculated rate constants. Appendix Proof of Corollary 1. Demonstration will be by construction of a Smith matrix with the specific property. Let X = (x&..~x,,) = X(0) be any Smith matrix and let p = (PI, ...,p,,) be an integer vector. Setpl = j for somej E 11, ..., K ) . Let c j k be the element in thejth row of X(0) having the smallest nonzero absolute value. Combine basis vector Xk with the others to obtain new matrix columns as follows

such that Ujk becomes the element of largest magnitude in thejth row. Find now the (different) value of k such that U j k has the smallest nonzero absolute value and repeat the process. Continue in this manner until, for some k, crj,#O,

uji=0, i = l ,

..., n,

i f k

(A.2)

Note that this simply amounts to the Euclidean algorithm for finding the greatest common divisor of the elements of the j t h row. Interchange the kthcolumn with the first and let the resulting matrix be called X(l). Now consider henceforth only columns 2, ...,n of X(l). Find a row j having a nonzero element in one of these columns, and set p2 = j . Apply the Euclidean algorithm as above, at the end interchanging the column containing the greatest common divisor with the second column, and call this matrix X(2). Continue recursively to obtain X(n). This matrix has the lower triangular submatrix L = (up,&xn. It now must be verified that Xcn)is a Smith matrix. X(")was created from the original X by performing only elementary column operations, hence X(n) is also a fundamental matrix solution of ( 6 ) . Adding multiples of columns to other columns does not change the value of the determinant of any of the n X n submatrices, and interchanging columns only changes the determinants in sign. Therefore, the operations applied do not alter the greatest common divisor of the determinants; hence the greatest common divisor of the determinants of the n X n submatrices of X(")is unity and X(n) is a Smith matrix. Proof of Corollary 2. Let X = (XI[ ...Ix,,) be a Smith matrix as in Corollary 1, and let L = ( h j k ) n x n be its lower triangular submatrix. For a l l j = 1, ..., n, such that A, < 0, set x, to -x, to

Computer Construction of Reaction Mechanisms

The Journal of Physical Chemistry, Vol. 97, No. 46, 1993 12105

give all positive elements on the diagonal of L. Since Ajj # 0, for k < j we may write by the division algorithm xjk

= qjk’jj + bjk

where q j k is an integer and 0 Ib j k xk

rendering

-

xj

equal to zero. If