Faraday Symp. Chem. SOC.,1984 19 137-147 MCSCF Gradient Calculation of Transition Structures in Organic Reactions BY FERNANDO AND ANDREA BERNARDI BOTTONI Istituto di Chimica Organica Universita di Bologna Bologna Italy AND JOSEPH AND MICHAEL J. w. MCDOUALL A. ROBB* Department of Chemistry Queen Elizabeth College London AND H. BERNHARD SCHLEGEL Department of Chemistry Wayne State University Detroit Michigan U.S.A. Received 16th July 1984 The applicability of MCSCF gradient methods to the calculation of transition structures and diradicaloid intermediates is discussed. It is shown how the diabatic surface model provides a useful criterion for the choice of the valence space in the MCSCF method and also provides useful qualitative information about the electronic rearrangement associated with various transition states.These ideas are then applied to the synchronous and asynchronous 1,3-dipolar cycloaddition of fulminic acid to acetylene. One of the major development areas in quantum chemistry in recent years has involved the computation and characterization of the intermediates and transition structures for model organic reactions. Often these structures are diradicaloid in nature and are not described even qualitatively at the SCF level. The investigation of this type of species has been facilitated by the development of the MCSCF method [for a comprehensive review see ref. (l)] and gradient-optimization techniques. In recent work we have been involved in the development of MCSCF2 gradient3 programs and in the subsequent analysis of the reaction profile using diabatic In the present communication we shall attempt to present a methodological strategy that can be used to characterize potential surfaces that may contain diradicaloid structures.The application of this strategy will be illustrated with a study of transition structures occurring in the 1,3-dipolar cycloaddition of fulminic acid to acetylene. The major difficulties in the application of MCSCF gradient methods to the computation of transition molecular structures lie (1) in the choice of configuration space of the MCSCF (i.e. the choice of the orbitals to be used in the CI expansion) and (2) in the determination of regions of the molecular potential surface in which to search for critical points.In this work we use a model for the description of saddle points based upon the surface of intersection of two diabatic surfaces (one associated with the reactants and one with products). Clearly the reference CI expansion must be chosen so that each diabatic surface is represented equally well. The diabatic surfaces can be computed approximately using the methods discussed in ref. (4) and (5). The subsequent geometry optimization is then performed using the methods of ref. (3). 137 CALCULATION OF TRANSITION STRUCTURES THEORETICAL METHODS The standard methods for the optimization of molecular geometries (e.g. the Newton-Raphson method with updating of the Hessian as used in the GAUSSIAN 80 series of programs) work quite well for diradicaloid transition structures (saddle points) provided one uses MCSCF gradient3 methods.However one must start the optimization procedure with a molecular structure where the force-constant matrix has one direction of negative curvature. Given the large number of molecular parameters in typical model organic reactions the location of these regions of the molecular potential surface can be quite difficult. Further it is obviously important not to make apriori assumptions about the nature of the reaction coordinate and thus miss possible transition structures. In recent work we have developed an MCSCF procedure4* for the computation of a diabatic surface in which the adiabatic surface of the reaction is obtained from the interaction of two diabatic surfaces.The main feature of this procedure is that each diabatic surface is associated with the specific bonding situation of either the reactants or products and the transition structure lies on the surface of intersection of the two diabatic surfaces. In ref. (5) we have shown that the geometries of the transition structure of the cyanate-isocyanate rearrangement the 1,2 sigmatropic shift in propene the S,2 reaction of H and CH and the addition of singlet methylene to ethylene correlate very accurately with the intersection of appropriate diabatic curves. In the case of the sigmatropic shift we were able to locate the transition structure a priori from preliminary diabatic-surface calculations. The significant feature of this approach is that not only does it furnish us with an excellent starting transition structure but also it provides some insight into the origin of the reaction barrier and conformational preference of possible reaction intermediates.The latter follows from the model used to define the diabatic surfaces themselve~.~?~ We now give some brief discussion of this model. One assumes that for a reaction X+Y the adiabatic surface can be described in terms of two sets of variables R and R,. We also assume that we can model this reaction with two diabatic surfaces E,(R, R,) and Ey(R, Ry).The surface Excorrectly describes the energy minimum of the species x and represents the bonding situation of the species x at all regions of configuration space (including the minimum of the surface E,).The definition of the two diabatic surfaces themselves involves (1) the use of fragment MCSCF orbitals of x and y and (2) a partition of the full valence CI space (CASSCF) into packets. Each packet consists of an isolated fragment configuration (either a Heitler-London or no-bond configuration in the language of valence-bond theory) plus all possible one-electron- transfer configurations. Thus for example in the cycloaddition of two ethylenes one packet (reactants) corresponds to two ethylene molecules in their ground states (a no-bond configuration) and the second (products) corresponds to two triplet excited Fthylenes coupled to an overall singlet (a Heitler-London configuration). The truncation of each packet at one-electron-transfer configurations will clearly begin to break down at small interfragment separations where more than one electron transfer will become very important and other locally excited configurations will begin to make large contributions.However the truncation at one-electron transfer has the feature that the packets are mutually exclusive and hence approximately orthogonal. With this definition at hand one can compute a diabatic surface4y5 by performing an MCSCF calculation on each packet separately. The possible transition structures correspond to minima on the surface of intersection of these diabatic surfaces. It is F. BERNARDI et al. apparent that the surface of intersection of two diabatic surfaces is essentially a conformational hypersurface thus the conformational minima (corresponding to various possible transition structures) allow a semi-classical analysis in terms of electrostatic polarization exchange repulsion and charge-transfer energies.In ref. (4) we have shown how this type of analysis can be performed for MCSCF wavefunctions. APPLICATIONS During the past two decades 1,3-dipolar cycloaddition reactions have become a general method for the synthesis offive-membered rings. However there is considerable controversy concerning the reaction mechanism.** There are two obvious alternatives a synchronous mechanism involving two-bond cycloaddition via an aromatic transition state or an asynchronous mechanism involving the formation of a diradicaloid intermediate. In the latter case one must find two transition states corresponding to the formation of the first and second bonds.While a synchronous cyclic transition structure can be determined by SCF methods,lo?l1 the possibility of an extended diradical intermediate can only be investigated using methods that include at least the two configurations necessary to describe the diradical. Recently,12 the addition of fulminic acid (1,3-dipole) to acetylene (dipolarophile) has been investigated using multireference CI methods. In these calculations a low-energy extended diradical intermediate was found. In the MCSCF calculations to be discussed here the synchronous and asynchronous processes can be treated with equal accuracy at the same level of theory. This problem is of considerable theoretical interest.First there are a very large number of possible 1,3-dipole/dipolarophile combinations with very different (computed) barriers for the synchronous process. l3 Secondly it should be apparent that the diabatic-surface method discussed above is capable of rationalizing the mechanistic controversy. In this formalism the preference synchronous/asynchronous is simply related to the relative positions of the two possible minima of the surface of intersection of the two diabatic surfaces corresponding to reactants and products. These facts in turn can be rationalized using the methods described in ref. (4). We shall begin our discussions with a qualitative discussion of the diabatic surface model as it applies to the 1,3-dipolar cycloaddition of fulminic acid to acetylene.Initially we are concerned with the choice of the valence orbital space (i.e.those active orbitals whose occupancy must be allowed to have values other than 2 or 0 and which will form the basis of the CI expansion used in the MCSCF calculations). This discussion is most conveniently given in terms of the orbitals of the isolated fragments and formulated in a valence-bond type of language. The fulminic acid can be thought of as having two allyl-like systems of three n orbitals referred to as the in-plane n (incipient 0 bonds) and out-of-plane n (incipient delocalized norbitals) sets as the two fragments approach each other. Similarly the acetylene has two sets of ethylene n orbitals which we shall refer to as the in-plane and out-of-plane set.The isolated fragments must have the orbital occupancy shown in fig. l(a) for the in-plane and out-of-plane n systems. For the product isoxazole the in-plane n-orbital occupancy must correspond to the promotion of one electron in each of the fragments from HOMO to LUMO. The unpaired electrons must be spin-coupled to a state of triplet spin within a fragment and these two triplet states subsequently spin-coupled to a singlet in order to describe two new 0 bonds as a valence bond Heitler-London configuration. For the out-of-plane n system of the product isoxazole the fragment n systems retain the configuration shown in fig. l(a) (ie.there is no requirement to uncouple and recouple the spins to describe the isoxazole n system). From the preceding argument it should be obvious that (in terms of the fragment CALCULATION OF TRANSITION STRUCTURES HCNO HCCH HCNO HCCH $3 -Fig.1. (a) No-bond and (b) Heitler-London configurations for the singlet-singlet and triplet-triplet diabatic surface for the addition of fulminic acid to acetylene. orbitals) one requires a valence space consisting of four in-plane n orbitals (the HOMO and LUMO of each fragment) in order to describe the product the transition state and the reactants with equal accuracy. Thus the valence space to be used in the MCSCF calculations should contain these four cr orbitals at least. Furthermore the product isoxazole the intermediate diradical and the transition state for the formation of the second bond should be dominated by the diabatic surface associated with the Heitler-London configuration (for the in-plane z orbitals) shown in fig.1 (b).We shall refer to this surface henceforth as the triplet-triplet surface. On the other hand the reactants (fulminic acid and acetylene) must be dominated by the diabatic surface associated with the no-bond configuration (for the in-plane z system) corresponding to the configuration shown in fig. 1 (a).We shall refer to this surface in subsequent discussions as the singlet-singlet diabatic surface. The transition structures for the synchronous formation of isoxazole and for the formation of the first cr bond in the asynchronous process should be associated with the intersection of these two diabatic surfaces. We can now outline our computational strategy for the characterization of the critical points on the surface for the cycloaddition of fulminic acid to acetylene.The orbitals for the initial MCSCF calculation for each geometry optimization were obtained from MCSCF calculations on the isolated fragments. These orbitals were then orthogonalized as follows (1) the core orbitals (doubly occupied in all reference CI configurations) the valence orbitals and the virtual orbitals were each symmetrically orthogonalized within each set and (2) the valence orbitals were then Schmidt orthogonalized to the core and the virtual orbitals subsequently Schmidt orthogonalized to the core and then to the valence orbitals. This procedure unambiguously defines the valence-orbital set. Note that the SCF orbitals are not suitable as starting orbitals because the virtual orbitals tend to be very diffuse with extended basis sets and the appropriate weakly occupied orbitals do not usually correspond to the lowest-energy occupied orbitals of an SCF calculation.All of our geometry optimizations are then carried out with a full CI in the space of the four in-plane 7z orbitals (corresponding to HOMO and LUMO of each fragment). In order to improve the energetics MCSCF calculations were also carried out (at the geometry just obtained) with a valence space that consisted of the HOMO-LUMO out-of-plane n orbitals as well corresponding to eight valence orbitals and 1764 configurations. This calculation should account for some of the dynamic correlation of the delocalized n system.The calculations were carried out at the STO-3G and 4-31G basis set level. At the STO-3G level each critical point was characterized by computing the hessian (by finite difference) in the subspace consisting of the interfragment geometrical parameters and F. BERNARDI et al. the CNO angle of fulminic acid which was strongly coupled to these variables. This hessian was then updated numerically in the 4-3 1G geometry optimizations. In fig. 2(a)-(f) we illustrate the geometries of the critical points located by the procedure described above. Table 1 contains the geometrical parameters of fulminic acid and acetylene obtained at the same level of theory. The absolute and relative energetics are summarized in table 2. It can be seen that there is a cyclic transition state [fig.2(b)] for the synchronous two-bond addition to give isoxazole [fig. 2(a)] an extended diradicaloid transition state [fig. 2(c)] for the formation of the first bond leading to a diradical intermediate which exists in a trans form [fig. 2(d)]or a cis form [fig. 2(e)] and a second diradicaloid transition state [fig. 201 which connects the cis form of the diradical intermediate [fig. 2(e)] with the product isoxazole. It can be seen that the STO-3G geometries are in quite good agreement with the 4-31G results indicating that the STO-3G basis is capable of giving a good qualitative description of the surface. The geometry of the cyclic synchronous transition state is in good agreement with that obtained in ref. (1 1) at the SCF level in the same basis.This is to be expected since this species is well described at the SCF level. The description of the diradicaloid region of the surface obtained in the present work is slightly different from that obtained in ref. (12) and deserves some comment. The geometry of the diradicaloid transition state [fig. 2(c)] that connects reactants and the diradical intermediate is in qualitative agreement with the structure obtained in ref. (12) at the UHF and 3x3 CI level (in the 4-31G basis). The major difference is that the interfragment C(3)-C( 1) distance is considerably shorter in the present compu- tations. However in the present calculations we find two different structures for the diradical intermediate a trans form [fig. 2(d)] and at slightly lower energy a cis form [fig.2(e)] that differ only in the arrangement of the C-C-H angles in the acetylene fragment.In the present work we have also searched for a cis form of the transition state shown in fig. 2(c). Since the calculation was clearly converging to the trans form it was abandoned Thus it seems likely that the transition state for the formation of the first bond connects the reactants with the cis form [fig. 2(e)] of the diradical intermediate and the trans form [fig. 2 (d)] represents a subsidiary minimum accessible uia inversion at C(2). The cis form of the diradical intermediate [fig. 2(e)] is again in good agreement with the structure found in ref. (12) where an improved estimate (obtained from large-scale multireference CI calculations) for the transition structure for the asynchronous process (formation of the first bond) was obtained by inter- polating geometries along a path connecting the cis form of the intermediate and the transition structure located at the SCF+ 3x3 CI level.An interpolation based on the trans form would have yielded a structure in better agreement with the present fully optimized transition state. Finally a transition structure has been obtained for the formation of the second bond in the asynchronous process [fig. 2(f)]. This lies energetically below the transition structure for the formation of the first bond in the asynchronous process. The total and relative (to reactants) energies are given in table 2 along with the CI results obtained from ref.(1 1) and (12) for comparison. It can be seen that the relative energetics for the synchronous process are in good agreement with those obtained from ref. (1 1) using a CI treatment that included replacements from two a’ and three aff orbitals. The MCSCF reaction barrier (synchronous process) of 26-28 kcal mol-l will be lowered by the inclusion of correlation in more a’ orbitals.11’l2 However it must be stressed that one must be careful in comparing CI calculations for the synchronous and diradicaloid regions of the surface. In the diradicaloid region of the surface SCF orbitals are far from optimum for use in the CI expansion even if a multireference expansion is used. Thus one must be certain that one is describing both regions of CALCULATION OF TRANSITION STRUCTURES Fig.2. Optimized geometries computed at MCSCF/STO-3G (no superscript) and MCSCF/4-3 1G (asterisk) (a) isoxazole (b) synchronous cyclic transition leading from fulminic acid and acetylene to isoxazole (c) asynchronous transition state corresponding to the formation of the first bond (i.e.connecting reactants and the trans diradical intermediate) (d) trans diradical intermediate (e)cis diradical intermediate and (f)asynchronous transition state connecting the cis diradical intermediate to the product isoxazole. 143 F. BERNARDI et al. 144 CALCULATION OF TRANSITION STRUCTURES Table 1. Optimized geometries computed at MCSCF/STO-3G and MCSCF/4-3 1G for fulminic acid (HCNO) and acetylene (HCCH) (atoms labelled as in fig.2) bond length/A molecule bond STO-3G 4-31G HCNO C(3)-H(3) 1.06 1.05 C(3I-W 1) 1.17 1.14 N( 1)-O( 11 1.29 1.26 HCCH C( 1)-H( 1) 1.07 1.05 C( 1)-C(2) 1.18 1.20 Table 2. Total and relative energies for the reaction of fulminic acid with acetylene MCSCFa MCSCFa MCSCF* structure STO-3G 4-31G 4-31G ref. (1 ref. (12)d total energies/hartreee HCCH +HCNO -241.3130 -244.1360 -244.2394 -244.5444 -244.6181 cyclic synchronous -241.3056 -244.0952 -244.1941 -244.5062 -244.5973 TSf [fig. Wl asynchronous first -241.2950 -244.0871 -244.1841 -244.603 1 TSf [fig. 2(c)l trans diradical -241.37 15 -244.1214 -244.2005 intermediate [fig. 2(d)] cis diradical -241.3746 -244.1246 -244.1977 -244.608 1 intermediate [fig. 2 (e)] asynchronous second -241.3636 -244.1 123 --TSf [fig.2cf)l isoxazole [fig. 2 (a)] -241.521 1 -244.2426 -244.3150 relative energies/kcal mol-1 HCCH +HCNO 0.0 0.0 0.0 0.0 0.0 cyclic synchronous 4.6 26.0 28.4 24.0 13.0 TSf [fig. 2(4l asynchronous first 11.3 30.7 34.7 -9.39 TSf [fig. 2(c)l trans diradical -36.6 9.1 24.4 intermediate [fig. 2(d)] cis diradical -38.6 7.2 26.2 -8.20 intermediate [fig. 2 (e)] asynchronous second -31.7 14.9 --TSf [fig. 2cf)l isoxazole [fig. 2 (a)] -130.6 -66.9 -47.4 -69.1 -63.4 a Full CI (CASSCF) in four a’ orbitals and four electrons. Full CI (CASSCF) in four a’ and four a” orbitals and eight electrons. CI using SCF orbitals all single and double excitations from two a’ and three a’’ orbitals. CI (two reference second-order Moller-Plesset) using SCF orbitals.1 hartree = 2625 kJ mol-l. f Transition state. F. BERNARDI et aZ. 145 the surface with equal accuracy. The MCSCF energetics clearly do not include any true dynamic correlation. Thus the exothermicity and the reaction barriers will not be accurately reproduced. Nevertheless the valence space has been chosen in this work so that the region of the synchronous transition structure and the region of the diradicaloid structures are described with equivalent accuracy. Thus one has some confidence in the relative energetics of the high-energy regions of the synchronous as against asynchronous pathways. From the MCSCF results presented in table 2 it can be seen that the transition state for the synchronous process is predicted to be lower than for the first step of the asynchronous process at all levels of computation.Furthermore the second transition state in the asynchronous processes lies at a lower energy than the first and the second barrier is predicted to be quite small. Thus both the synchronous and asynchronous processes could be consistent with retention of stereochemical information (in the case of substituted adducts given the same relative orientation of approach). Finally we turn to an a posteriori rationalization of some of these results in terms of the diabatic-surface model. This is instructive because it is indicative of the success of the model when used to locate transition structures apriori and it also complements the qualitative discussion given previously.In fig. 3-5 we give cross-sections through the diabatic surfaces chosen to pass through the cyclic synchronous transition structure (fig. 3) the first asynchronous structure (fig. 4) and the trans diradical -240.80 -241.00 -24140 . 2.0 2.5 3.0 3.5 4.0 RXdB Fig. 3. Diabatic surface cross-section passing through the synchronous cyclic transition state (-) singlet-singlet (--) triplet-triplet. intermediate minimum (fig. 5). The figure in each case shows the diabatic energies obtained along the coordinate (RXN)corresponding to a rigid dissociation along the line connecting N(l) of the fulminic acid fragment and the centre of the C(l)-C(2) bond. In fig. 3 (cyclic synchronous transition state) one observes that the singlet-singlet diabatic curve is very flat and begins to become repulsive only in the region of the transition state (RxN= 2.74 A) itself.In contrast in fig. 4 (corresponding to the first diradicaloid transition structure) the singlet-singlet diabatic surface quickly becomes repulsive and intersects the attractive triplet-triplet surface in the region of the transition structure (RxN = 3.41 A). The intersection of the singlet-singlet and triplet- triplet surfaces in the synchronous process occurs in a region where the truncation at one-electron transfer is no longer valid. However in fig. 5 (diradical intermediate) CALCULATION OF TRANSITION STRUCTURES -24 0.80 1 t -241.00 irl -241.20 -241.40 t 2.0 2.5 3.0 3.5 4 .O RXdB Fig.4. Diabatic surface cross-section passing through the asynchronous transition state for the formation of the first bond (-) singlet-singlet (--) triplet-triplet. -240.60--240.80 -a 2 -241.00 s -. t 4; 9 I -241.20 : -241.40 Y 2.0 2.5 3.0 3.5 4.0 RXdB Fig. 5. Diabatic surface cross-section passing through the trans diradical intermediate (--) singlet-singlet (--) triplet-triplet. the geometry of the fragments is now such that the triplet-triplet diabatic surface lies below than the singlet-singlet for all values of the interfragment separation. Thus the synchronous transition state and the first asynchronous transition state correspond to regions of intersection of the diabatic surfaces. However the second asynchronous transition state (from the diradical intermediate to isoxazole) lies (like the diradical intermediate itself) on the triplet-triplet diabatic surface.The synchronous two-bond cycloaddition appears to be feasible in this situation because of the fact that the singlet-singlet diabatic surface is not strongly repulsive due to the stabilization of the extended n system. F. BERNARDI et al. CONCLUSIONS It is clear that SCF methods cannot provide an adequate description of both the synchronous and asynchronous pathways for 1,3-dipolar cycloaddition reactions whereas MCSCF methods will provide an accurate description of both processes. For the asynchronous process MCSCF is necessary both to describe the avoided intersection of the singlet-singlet and triplet-triplet surfaces and because the intermediate is diradicaloid.The use of the diabatic surface model provides a useful scheme on which to base one’s selection of the valence space and also provides some a posteriori rationalization of the topology of the surface. Note added in proof We have also located (at the STO-3G level) the transition state connecting the cis [fig. 2(e)]and trans [fig 2(f)]diradical intermediates. The energy of this transition state is 7.5 kcal mol-l above that of the trans intermediate. All computations were run on the CDC 7600 or the CRAY 1s computers at the University of London national computing centre and on the VAX 11/750 at Queen Elizabeth College. The cooperation of both computer centres is gratefully acknowledged.The MCSCF gradient programs have been installed as part of the GAUSSIAN 8014suite of codes. We also acknowledge the financial support of NATO under grant no. RG 096.81. A. B. thanks the Royal Society for the award of a visiting fellowship. J. Olsen D. L. Yeager and P. Jorgensen Adv. Chem. Phys. 1983 55 1. R. H. A. Eade and M. A. Robb Chem. Phys. Lett. 1981,83 362. H. B. Schlegel and M. A. Robb Chem. Phys. Lett. 1982 93,43. F. Bernardi and M. A. Robb Mol. Phys. 1983,48 1345. F. Bernardi and M. A. Robb J. Am. Chem. Soc. 1984 105 54. F. Bernardi M. A. Robb H. B. Schlegel and G. Tonachini J. Am. Chem. SOC.,1984 106 1198. F. Bernardi A. Bottoni and M. A. Robb Theor. Chim. Acta 1984 64 259. * R. Huisgen Angew. Chem. 1963 2 565. R. A.Firestone J. Org. Chem. 1968 33. 2285. lo D. Poppinger Aust. J. Chem. 1976 29 465. l1 A. Komornicki J. D. Goddard and H. F. Schaefer J. Am. Chem. SOC.,1980 102 1763. l2 P. C. Hiberty G. Ohanessian and H. B. Schlegel J. Am. Chem. SOC.,1983 105 719. l3 M. Sana and G. Leroy J. Mol. Struct. 1982 89 147. l3 GAUSSIAN 80 J. S. Binkley R. A. Whiteside R. Krishnan R. Seeger D. J. DeFrees H. B. Schlegel S. Topiol L. R. Kahn and J. A. Pople Quantum Chemistry Program Exchange 1981 13 406.