Cundari Th.R. -- Computational Organometallic Chemistry-0824704789

158

Pierloot

47.AW Ehlers, S Dapprich, SF Vyboishchikov, G Frenking. Organometallics 15:105– 117, 1996.

48.BJ Persson, BO Roos, K Pierloot. J Chem Phys 101:6810–6821, 1994.

49.K Pierloot, E Tsokos, LG Vanquickenborne. J Phys Chem 100:16545–16550, 1996.

50.BO Roos, K Andersson, MP Fu¨lscher, L Serrano-Andre´s, K Pierloot, M Mercha´n, V Molina. J Mol Struct (Theochem) 388:257–276, 1996.

51.K Pierloot, E Van Praet, LG Vanquickenborne, BO Roos. J Phys Chem 97:12220– 12228, 1993.

52.A Ceulemans, GA Heylen, LF Chibotaru, TL Maes, K Pierloot, C Ribbing, LG Vanquickenborne. Inorg Chim Acta 251:15–27, 1996.

53.A Ceulemans, LF Chibotaru, GA Heylen, K Pierloot, LG Vanquickenborne. Chem Rev 100:787–806, 2000.

54.C De Graaf, F Illas, R Broer, WC Nieuwpoort. J Chem Phys 106:3287–3291, 1997.

55.C De Graaf, R Broer, WC Nieuwpoort. Chem Phys Lett 271:372–376, 1997.

56.C De Graaf, IDR Moreira, F Illas, RL Martin. Phys Rev B Cond Mat 60:3457– 3464, 1999.

57.M Flock, K Pierloot. J Phys Chem A 103:95–102, 1999.

58.J Sugar, J Corliss. J Phys Chem Ref Data 8:1–62, 1979.

59.J Sugar, J Corliss. J Phys Chem Ref Data 10:1097–1174, 1981.

60.CE Moore. Atomic Energy Levels. Vol. III. Washington, DC: Circular of the National Bureau of Standards 467, 1958.

6

Quantitative Consideration of Steric Effects

Through Hybrid Quantum Mechanics/

Molecular Mechanics Methods

Feliu Maseras

Universitat Auto`noma de Barcelona, Barcelona, Catalonia, Spain

1. INTRODUCTION

Steric effects have been largely absent from the spectacular progress experienced by computational organometallic chemistry in the last decades. There is good reason for this. The methodological and computational struggle to properly describe the properties at transition metal centers leaves little space for the introduction of the bulky ligands responsible for steric effects. This situation is, however, currently changing, in part because of the entry of hybrid quantum mechanics/ molecular mechanics (QM/MM) methods into this field of chemistry.

The power of QM/MM methods is based directly on that of pure QM and pure MM methods. There are standard QM methods that describe reliably all chemical features of small molecular systems. There are standard MM methods that describe reliably some chemical features of large molecular systems. Both types of methods, pure QM [either Hartree–Fock (HF)–based or density functional theory (DFT)] and pure MM, are steadily expanding the range of systems to which they can be applied, and up-to-date notice of this can be found in other chapters of this same volume. The QM/MM approach is simpler: The chemical

159

160

Maseras

system is divided in two regions, and the more convenient method is applied to each region. In this way, the heavy computational cost of QM methods can be concentrated only in the regions where it is strictly required while keeping an overall correct MM description for the rest of the system. The potential application of QM/MM methods to organometallic chemistry is enormous, because the electronic complexity of the systems is usually concentrated in a small region, namely, the metal and its immediate environment.

This chapter is intended for an audience of computational organometallic chemists interested in the practical use of hybrid QM/MM methods. Because of this the description of the methodological details will be kept to a minimum, condensed in the second section. Similarly, the chapter is not intended to be a review of published applications, which can be found in another recent review

(1). Instead, this chapter illustrates a series of practical aspects of QM/MM calculations that make them different from other, traditional QM or MM approaches. These features are shown mostly through the presentation of selected aspects of different examples of calculations, some of them carried out specifically for this text, some of them taken from previous publications, mostly by the author.

Most of the applications presented use the integrated molecular orbital/ molecular mechanics (IMOMM) method (2), which therefore will be briefly described in the next section. Afterwards, two sections will discuss particular aspects of the calculation setup, and another section will present specific ways to analyze QM/MM results. A final section will offer concluding remarks.

2.INTEGRATED MOLECULAR ORBITAL/MOLECULAR MECHANICS METHOD WITHIN THE CONTEXT OF QUANTUM MECHANICS/MOLECULAR MECHANICS METHODS

Hybrid QM/MM methods already have a certain history of their own in computational chemistry (3–8). Early work in this field has been on the introduction of solvation effects, with special focus on biochemical systems. The general application of this approach to transition metal systems, where there are often chemical bonds across the frontier between QM and MM regions, has been more recent. The IMOMM method, proposed in 1995 (2), has been remarkably successful, as shown by the large number of applications, concentrated mostly in transition metal chemistry (9–28), although not exclusively (29–31). The method has also been the starting point of other QM/MM methodological developments (32–36). One must cite in this regard the IMOMO and ONIOM methods. The IMOMO method (32) is the extension of the method to the use of two different-quality QM descriptions. The ONIOM method (33,34) is essentially a generalization that encompasses both the IMOMM and IMOMO methods, with the significant addition of the possibility of using more than two layers. The reason why the

Quantifying Steric Effects via QM/MM Methods

161

label IMOMM (instead of the more general ONIOM) is used throughout this chapter is fundamentally practical: the program used in most of the calculations presented was that of the original IMOMM implementation.

Good general discussions on hybrid QM/MM methods can be found in presentations of methodological novelties (8,34,37) and also in recent reviews of these methods (1,38–41). Because of this, the discussion here will be very brief.

The main difference between the current implementations of IMOMM, IMOMO, and ONIOM and the majority of other available QM/MM methods is related to the handling of the interaction between the QM and the MM regions. In principle, in any hybrid QM/MM method the total energy of the whole system can in all generality be expressed as:

Etot(QM, MM) EQM(QM) EMM(MM) Einteraction(QM/MM)

where the labels in the subscript refer to the type of calculation and the labels in parentheses correspond to the region under study. Both the QM and MM methods can in principle compute the interaction energy between the QM and MM regions, and the previous expression becomes:

Etot(QM, MM) EQM(QM) EMM(MM)

EQM(QM/MM) EMM(QM/MM)

The energy expression in a general hybrid QM/MM method thus has four components. Two of them correspond simply to the pure QM and MM calculations of the corresponding regions. And the other two correspond to the evaluation of the interaction between both regions, in principle at both computational levels. Different computational schemes are defined by the choice of the method to compute the EQM(QM/MM) and EMM(QM/MM) terms.

One of the defining characteristics of current implementations of IMOMM and derived methods is the neglect of the EQM(QM/MM) term. The neglect of this term obviously introduces an error in the reproduction of experimental reality. But it simplifies the calculation enormously, from a technical point of view, and it leads to an easier interpretation of the results.

This EQM(QM/MM) term is usually critical in solvation problems, because one of the points of interest is precisely how the quantum mechanical properties of the solute are modified by the presence of the solvent. In the case of a transition metal complex, this term would account mainly for the electronic effects of the ligand substituents on the metal center. A common way to introduce the EQM(QM/ MM) term is to put electrostatic charges in the positions occupied by the MM atoms, introducing in practice a term in the monoelectronic Hamiltonian (8). One problem with this kind of approach is the choice of the electrostatic charges, which is by no means trivial, since its validity is usually confined to the consis-

162

Maseras

tency of a force field. A more serious problem appears when there are chemical bonds between the QM and MM regions, where this approach breaks down in the proximity of the interface and needs to be reformulated. A more elaborate scheme to answer this problem has been proposed through the introduction of localized orbitals, but so far it has been applied mostly to organic systems (42,43). As mentioned earlier, IMOMM neglects this EQM(QM/MM).

The EMM(QM/MM) term, which is considered explicitly by IMOMM and derived methods, accounts for the direct effect of the atoms in the QM region on the MM energy of the system. This term is usually critical in transition metal complexes, because it accounts for the geometrical constraints introduced by the presence of the metal center in the arrangement of the ligands. In other words, it is related mostly to the steric effect of the ligand substituents. The EMM(QM/ MM) term is usually introduced through the parameterization of the QM atoms with the same force field used in the MM region. This parameterization is much simpler than it would be in the case of a pure MM calculation, because the MM part of the QM/MM calculation neglects the interactions within the QM region.

The presence of chemical bonds between the QM and MM regions poses a problem to the use of hybrid QM/MM methods, with different approaches taking different solutions. The particular method applied by IMOMM and derived methods, the introduction of additional link atoms to saturate the dangling bonds, will be discussed in more detail in Section 3.4 of this chapter.

It must also be mentioned that IMOMM and derived methods involve a full multistep optimization (1). This means that the geometry of both the QM regions and the MM regions is modified to minimize the total energy, with the result that the final geometry corresponds neither to the optimal QM arrangement nor to the optimal MM arrangement. Simpler one-step methods, where the MM geometry is optimized on a frozen QM geometry, also have their value but will not be discussed here.

The IMOMM method can in principle be applied to any combination of QM and MM methods. Both the QM and the MM level of the calculation are indicated in this text through a compact terminology of the type IMOMM(QM level:MM level).

3. DEFINITION OF THE COMPUTATIONAL LEVEL

The chemist interested in performing a hybrid QM/MM calculation must make a number of choices that may affect substantially the result of the computation. These choices, ranging from the obvious to the subtle, will be discussed in this section.

3.1.Quantum Mechanics/Molecular Mechanics Partition

A fundamental decision to be taken in the planning of a QM/MM calculation is that of the QM/MM partition: Which atoms are going to be included in the QM

Quantifying Steric Effects via QM/MM Methods

163

region and which atoms are going to be included in the MM region? The guiding idea for this choice must be to use a QM region as small as possible yet containing all the fundamental interactions that cannot be described by the MM method. Chemical knowledge is the main guiding line for this choice, although calibration tests will be necessary in doubtful cases. These tests should ideally be carried out through comparison of the hybrid QM/MM results with those of pure QM calculations for the whole system. This would provide a more reliable criterion than comparison with experiment, because there are a number of reasons unrelated to the QM/MM partition why QM/MM results may differ from experimental data, like inaccuracies in the QM description or the presence of solvent or packing effects in the experimental data. Comparison with experiment is nevertheless still a useful criterion, especially when the agreement is satisfactory.

An interesting example showing the importance of the QM/MM partition can be found in a study containing IMOMM(MP2:MM3) calculations on potentially agostic Ir(H)2(PR3)3 complexes (20). An agostic interaction is the intramolecular interaction that takes place within one complex between the metal center and a C–H bond of one the ligands (44). These Ir(H)2(PR3)3 complexes have an empty coordination site at the metal and CH bonds in the ligands, but present agostic interactions only for certain combinations of R substituents in the phosphines.

In particular, we will discuss here the results on Ir(H)2(PCy2Ph)3 . In this

complex, X-ray data (20) indicate the existence of an agostic distortion, with one of the Ir–P–C angles being 100.9°, significantly different from the standard

tetrahedral angle of ca. 109°. In the IMOMM(MP2:MM3) calculation of this species, two different partitioning schemes between the QM and MM domains were applied (Fig. 1). In model I, all atoms not directly bound to the metal were

FIGURE 1 Different QM/MM partitions used in the IMOMM calculations of Ir

(H)2(PCy2Ph)3 . Atoms in the QM part are shown in black.

164

Maseras

calculated at the MM level, the QM part therefore being Ir(H)2(PH3)3 . In model II, the QM part of the phosphine included the agostically distorted chain, the QM part thus being Ir(H)2(P(Et)H2)3 . The essential difference between both models is thus in the description of the potentially agostic C–H bond.

The geometry of Ir(H)2(PCy2Ph)3 was fully optimized at the IMOMM(MP2:MM3) level with both QM/MM partitions I and II. The

IMOMM(MP2/I:MM3) calculation gave results in qualitative agreement with X- ray data, with the largest discrepancy being in the value of 105.6° for the agostic Ir–P–C angle. This was larger than the experimental X-ray value of 100.9° but

already smaller than the computed average Ir–P–C value at this phosphorus center of 109.9°. This result is interesting because it proves that the bulk of the

ligands alone is able to push one of the C–H bonds of the cyclohexyl group to the proximity of the metal, but it is still somehow removed by 4.5° from the experimental value. Use of the more elaborate model IMOMM(MP2/II:MM3),

with the C–H bond in the QM part, led to results much closer to the X-ray values. The Ir–P–C bond angle improved to 99.8°, only 1.1° from the experimental

value. This change proves that the interaction between the Ir center and the C– H bond is not properly reproduced by the MM3 force field, and must therefore be included in the QM region. The problem is obviously that the force field is not parameterized to describe agostic interactions.

The definition of the QM and MM regions is therefore critical for the validity of the QM/MM calculation. The smaller the QM region, the more affordable the computational cost, but care must be taken not to leave any critical electronic interactions out of the QM region. The importance of the choice of the QM/MM partition must, however, not hide the fact that both the QM and the MM descriptions must describe with sufficient accuracy interactions within the respective regions.

3.2.Quantum Mechanics Level

There is a wide variety of QM levels that can be applied to organometallic compounds, ranging from semiempirical methods like PM3(tm) to high-level methods such as multireference-configuration interaction, passing through all the derivations of HF and DFT methods. Information about these methods can be found in other chapters of this book. The choice of the most appropriate QM method for each chemical problem is by no means trivial. One of the main features of the use of hybrid QM/MM methods is the use of a small QM region, allowing the use of QM levels that would be unaffordable with pure QM calculations for all the system. The QM level chosen must nevertheless be sufficiently accurate to describe all significant interactions within the QM region, or else the whole QM/MM calculation will fail. In what follows, we will discuss one example showing how this can be critical.

Quantifying Steric Effects via QM/MM Methods

165

This example concerns calculations on heme species. Understanding the features of this type of complexes is important because heme groups are at the active center of a number of biochemically very relevant proteins and enzymes, like hemoglobin, myoglobin, cytochromes, peroxidases, and catalases (45). However, its theoretical study has been hindered so far by the large size of the system, containing a porphyrin ring (4 nitrogen atoms plus 20 carbon atoms) attached to an iron center. Because of this, it would be appealing to apply a QM/MM method to reduce the computational effort. This requires introducing the QM/MM partition within the porphyrin ring, which will obviously worsen the modeling of the ring aromaticity. The question of the magnitude of the effect of the introduction of the QM/MM partition on the description of the electronic properties of iron has no straightforward answer, and required performing a series of systematic test calculations (22). The conclusion of that study was that the partition presented in Fig. 2, cutting the QM part of the heme group to [Fe(NH(CH)3NH)2], leads to satisfactory results. This study nevertheless had an interesting spinoff concerning the validity of the RHF description that merits discussion here.

The particular complex discussed here is [Fe(P)(Im)(O2)] (P porphyrin, Im imidazole) (22). This species is a model for the active center of oxygen transport proteins hemoglobin and myoglobin, and both in the biological and bioinorganic systems is able to bind oxygen. In spite of that, IMOMM(RHF: MM3) calculations on this system did not yield any stable minimum with the

˚

oxygen bound to the iron, the optimal Fe–O distance being above 3.5 A. A

˚

weakly bound state (with bond distances around 2.2 A) could be obtained by

FIGURE 2 QM/MM partition used in IMOMM calculations of systems involving the heme group. Atoms in the QM part are shown in black.

166 Maseras

using a minimal basis set on oxygen, but this value is far from the experimental

˚

1.746 A found in related bioinorganic models and can almost certainly be attributed to a computational artifact associated with basis set superposition error. Therefore, in view of these results, one can see that the calculations fail to reproduce, even at a qualitative level, experimental reality. The QM/MM partition was proposed to be responsible for this failure.

However, this hypothesis was proved wrong by additional calculations. When IMOMM(Becke3LYP:MM3) calculations were carried out with the same

˚ ˚

QM/MM partition, the Fe–O distance became 1.759 A, only 0.01 A away from the X-ray value. Therefore, the inaccuracy of the IMOMM(RHF:MM3) calculation had nothing to do with the fact that part of the system was described with an MM method, but was due to the failure of the RHF method in the description of the iron–oxygen interaction in this particular system.

3.3.Molecular Mechanics Level

As with the QM level, the MM level utilized can also affect decisively the outcome of the calculation. However, the choice of an appropriate MM level has some complications that were absent in the case of the choice of the QM level described earlier. The first of these complications is that there is no clear hierarchy of MM methods. In contrast with QM methods, all MM methods have similar computational costs, and the differences between them are concentrated mostly in the type of system for the which they are parameterized.

An additional problem, of a technical nature, is that while QM methods are usually available as options of a single program, MM force fields are usually available only in independent programs. On one hand, this poses a serious limitation in the comparison of their performances, although some efforts have been reported (46). On the other hand, this requires the programmer interested in building a QM/MM code to make a specific interface for each of the force fields, and to have access to the source code of each of them. An unfortunate outcome of this situation is that most applications of QM/MM methods to organometallic chemistry have been carried out with a single force field, namely, MM3 (47). This force field, devised especially for organic systems, seems appropriate for the description of steric effects, and in fact the comparison with experiment of the results obtained in its application in QM/MM methods is mostly correct. However, it would be highly desirable to carry out systematic comparisons with the performance of other force fields, and so far this has seldom been possible. The most significant results in this direction are probably the two independent sets of QM/MM calculations on the particular problem of olefin polymerization via homogeneous catalysis. One set of calculations (25,26) was carried out with the MM3 force field; the other set (27,28) was carried out with the AMBER force field. The results were in remarkable agreement, indicating that, at least in this particular case, the choice of one of the two force fields was not critical.

Quantifying Steric Effects via QM/MM Methods

167

A second possible approach to the topic of the quality of the MM description is the tuning of some parameter of the force field to the particular problem under study. Although this has seldom been done in applications of the IMOMM method to transition metal chemistry, because it can obscure interpretation of the result, there is at least one particular example where it has been helpful (12).

The problem arose from the fact that IMOMM(Becke3LYP:MM3) calculations of a series of compounds containing chlorine ligands yielded abnormally poor results. Among these, we are going to discuss here the case of the Ir(H)2Cl(P-

(tBu2Ph))2 complex (Fig. 3). This complex has an experimental asymmetry (48) between the two H–Ir–Cl bond angles (156° and 131°), which was not repro-

duced by the IMOMM(Becke3LYP:MM3) calculation with the standard MM3 parameters (bond angles of 147° and 146°). This problem was solved from the

observation that the van der Waals radius of chlorine used in the MM3 force field, which is fundamental in defining its steric activity, is defined for organic systems, where chlorine must have less anionic character than in inorganic compounds.

The relationship between the van der Waals radii of organic and inorganic

chlorine was explored through high-level calculations on the CH3–Cl He and Na–Cl He model systems, and it was found that the van der Waals radius for

˚

inorganic chlorine is larger by 0.4A. This increase was introduced in the corresponding MM3 radius, and the resulting IMOMM(Becke3LYP:MM3) calculations provided results in much better agreement with experiment. In particular,

for the iridium complex described earlier the computed H–Ir–Cl bond angles became clearly different, with values of 162.6° and 122.0°.

3.4. Bond Distances of Connecting Atoms

One of the subjects that has taken a good deal of space in methodological discussions of the design of QM/MM methods is the way to deal with the connection between the QM and MM regions when there are chemical bonds across the boundary. In a number of methods, including IMOMM and derived methods, additional link atoms (usually hydrogen) are introduced in the QM calculation

FIGURE 3 Two views of the geometry of the Ir(H)2Cl(PtBu2Ph)2 complex.