On the Importance of Considering Multinuclear Metal Sites in Homogeneous Catalysis Modeling

In this short review, we provide an account of a number of computational studies of catalytic reaction mechanisms carried out in our groups. We focus in particular on studies in which we came to realize during the course of the investigation that the active catalytic species was a bimetallic complex, rather a monometallic one as previously assumed. In some cases, this realization was in part prompted by experimental observations, but careful exploration based on computation of the speciation of the metal precursor also provided a powerful guide: it is often possible to predict that bimetallic species (intermediates or transition states) lie lower in free energy than a priori competitive monometallic species. In this sense, we argue that in organometallic catalysis, the rule whereby “two is better than one” turns out to be relevant much more often than one might expect.


Introduction
Quantum chemistry is nowadays one of the most powerful tools for studying reaction mechanisms in homogeneous catalysis. Improved electronic structure methods and ever faster implementations thereof, together with computers that are more and more powerful mean that the range of problems that can be tackled continues to expand, while at the same time, more accurate results can be obtained for a given level of complexity of the system. As a result, computation is also becoming more and more widely used in mechanistic studies. These developments, as well as many outstanding applications, have been discussed in a number of recent reviews [1][2][3][4][5][6][7][8]. While this ever-larger role for quantum chemistry is welcome in many respects, it should also be borne in mind that such computational studies tend to have their own challenges that must be carefully considered when assessing the results.
In this account, we wish to address one issue that has arisen a number of times in studies from our two groups, and which we feel is important enough to warrant highlighting as a topic in its own right. We wish to emphasize already here that this is probably not new knowledge to the experienced practitioners. Nevertheless, even though this question is by now well-known to us, we still tend to overlook it on occasion, and believe that this may happen to others also. In that sense, we hope that this essay can serve as a useful reminder of the common nature of this issue.
The problem is one example of the more general problem relating to understanding catalytic reaction mechanisms whereby the active catalyst and/or the most abundant form of the catalyst in the reaction medium are not necessarily the same as the catalyst precursor initially introduced in the reaction vessel. The 'Halpern paradox' is a well-known phenomenon [9,10] in catalysis, whereby the most stable form of the catalyst under reaction conditions does not necessarily contribute to the main turnover cycle (and indeed, in some cases, is likely not to, hence the name 'paradox'). More generally, it quite often occurs that the most stable form of the catalyst differs from the pre-catalyst due to some pre-transformation. The specific topic covered here is an example of the latter case: although many catalysts 1 3 are introduced as monometallic precursors, a non-negligible number of catalytic reactions seem to involve dimeric metal complexes as active species. We do not claim that this will always be so (so the principle reviewed here is in no ways a paradox), but simply that it will be the case much more often than one might expect.
When performing a computational study of a metal-catalyzed transformation, the first step is to establish the nature of the active catalyst. In the early days of computational studies of catalysis, this information was usually taken 'on trust' as information coming from the experimental collaborator, or as general knowledge. As computational studies have become more and more ambitious, though, the active catalyst in the reaction is nowadays in many cases not established, and the computational mechanistic investigation must therefore either start with an assumption regarding the nature of this species, or set out to establish it from first principles. This latter procedure can be described as being the "computational speciation" study of the catalyst: one suggests possible transformations that the pre-catalyst may undergo in presence of solvent, reagents, additives or impurities, and sets out to compute the relative standard free energy of these species, to work out which are the most stable. As always in computational studies, while this procedure can in principle be predictive, the existence of significant error bars on free energy predictions means that it is wise to take into account as much experimental evidence as possible to establish whether the computational suggestion is compatible with known facts.
We have noticed from examples in our own work that in the absence of direct experimental evidence to the contrary, in cases where the pre-catalyst is a monometallic complex, the default assumption for the active species is usually that it too is monometallic. Even when computational speciation is performed, it tends to limit itself to monometallic species. This is perhaps not too surprising: textbook reaction mechanisms for elementary steps in organometallic chemistry involve single metal centers placed as the 'conductor' at the center of the catalytic orchestra, and many catalyst precursors and isolated intermediates involve single metals. Nevertheless, as is also well known, the free energy difference between monometallic complexes and corresponding multimetallic clusters is quite often small, so one should expect that depending on conditions, the most stable species might involve multiple centers playing together like in chamber music. Such polymetallic active species are trivially identified in the cases where the pre-catalyst contains multiple metals, as for example in the case of dirhodium catalysis [11], but are much less easy to predict when this is not the case.
As discussed above, we are aware that other researchers might have similar experiences and that there have been discussions on mono-versus bimetallic mechanisms in recent years in, for example, gold and palladium catalysis [12,13]. Nevertheless, in our overview below, we focus mainly on cases taken from our own research, because in these cases, we have the advantage that we know about the "order of discovery" of the work (rather than just about the "order of exposition" in the published paper), so we also know that most of these cases, the "bimetallic" mechanism was a late hypothesis that occurred only after a frustrating amount of time spent trying to account for observations based on monometallic mechanisms. Also, we will not put much effort on describing the experimental backgrounds of the discussed reactions nor their importance and synthetic value; instead, the focus will be on the modeling aspects.
arising from modes with extremely low frequencies, or from modes erroneously predicted to be imaginary. Use of the quasi-harmonic approximation can be helpful in these cases [14,15].
A specific comment on the role of entropy in dimer formation is perhaps useful. At first sight, when considering a reaction step 2A → A 2 , one might expect that this would be entropically disfavored as described by the simple 'rule of 10' (whereby ΔG = ΔH − TΔS is roughly ten kcal/mol less favorable at room temperature than ΔH due to the entropic cost associated with loss of translational and rotational degrees of freedom upon combining two species to form a single one). It should however be pointed out that for metal complexes, this expectation is too simplistic and the entropic contribution to a dimer-forming step can be quite different than that outlined here. For example, dimer formation can also involve release of a ligand, such that it is better written as 2(ML n ) → (ML n-1 ) 2 + 2L. In this case, dimerization will actually be favorable in entropic terms, as two molecules are combined to form three. All of our calculations in the studies mentioned below take entropy into account in the most accurate way possible so such effects are handled (though the detail will not be discussed below).

Pd-Catalyzed Alkene Isomerization
An early example of the problem with mono-vs. bimetallic catalyst arose when studying cis-trans isomerization of alkenes by Pd(II) [16]. The used catalyst (Scheme 1) is monometallic and no obvious need for a second metal arose from initial consideration of possible mechanisms. Indeed, a reaction path for isomerization after complexation to a single Pd(II) center was rapidly found in the computational study, though with a free energy barrier that was too high compared to the observed reactivity.
The initial temptation was to blame this on computational error or on missing features of the computational model other than its monometallic character, but changes in level of theory, attempts to include co-reactants or changes in the ligands around the metal did not yield any lowering of the barrier. Then the experimental kinetics provided a hint: the order of the reaction with respect to the quantity of catalyst used was not exactly 1, which led to the idea that a bimetallic complex might be involved. After some work, a revised (and unexpected) mechanism involving nucleophilic attack on an alkene coordinated to one metal by chloride connected to the other metal was found (Scheme 1). This model immediately led to a predicted free energy profile that was in good agreement with the experimentally derived one, even with DFT methods. It should be noted that the dimerization of PdCl 2 L 2 to form Pd 2 Cl 4 L 2 leads to release of ligand L, so is an example of the type of process 2(ML n ) → (ML n-1 ) 2 + 2L mentioned in the previous section. The calculated ΔG agreed quite well with the one inferred from the experimental kinetics) [16], suggesting that our computational protocol was able to handle the entropic aspects of speciation with sufficient accuracy.

Zr-Catalyzed Amide Formation
The second example we introduce here is the case of the direct condensation of carboxylic acids and primary amines (Scheme 2) [17]. Despite the importance of the direct catalytic amidation reaction in synthetic organic chemistry, the detailed mechanisms of this type of reactions were not well understood. A DFT mechanistic

+ HOH
study was performed to elucidate the mechanism of the ZrCl 4 -catalyzed reaction, using phenylacetic acid and benzyl amine (R 1 = R 2 = Bn) as model substrates.
Kinetics studies and NMR spectroscopy, in combination with precedence from the literature, provided evidence that the active complex contains up to two carboxylate ligands per Zr ion, coordinated bidentately to the metals. Moreover, the experiments suggested that the active catalyst complex contains one amine ligand per Zr ion. With these boundary conditions, a computational speciation study was conducted first showing that the dinuclear complexes are much more stable than mononuclear ones. Complex 4 was found to be as much as 42 kcal/mol lower in energy than complex 5 (Fig. 1). Similar dinuclear Zr complexes of the type [Zr(OR) 4 (RNH 2 ) 2 ], with two bridging alkoxide ligands, have precedents in the literature [18][19][20].
From complex 4, the catalytic mechanism proposed on the basis of the calculations is shown in Fig. 2. First, a nucleophilic attack of BnNH 2 at the carboxylate carbon gives the tetrahedral intermediate Int1 with a barrier of 10.7 kcal/mol (TS1). Next, an external amine enters the catalytic cycle and abstracts a proton from Int1 to form Int2, which undergoes a C − O bond cleavage and a proton transfer (TS2), leading to the generation of Int3. TS2 is 21.6 kcal/mol higher in energy than 4 and constitutes the rate-determining step of the catalytic cycle. Protonation of the hydroxyl group in Int3 by an external phenylacetic acid releases a water molecule and the amide product while regenerating catalyst complex 4.
Very interestingly, the reaction starting from mononuclear complex 5 resulted in the same step sequence as the dinuclear complex, but with a much lower barrier, only ca 10 kcal/mol, which is not compatible with the experimental conditions. The fact that the bimetallic complex 4 is more than 40 kcal/mol more stable than the monometallic complex 5 shows that the latter will not be present in the solution and cannot contribute to the catalysis, despite the low barrier to reaction that would result were this species to contribute to the mechanism. This example demonstrates thus that if one is not careful in the computational speciation study, and one thereby fails to exclude complex 5 on the basis of its high energy relative to 4, wrong conclusions can be drawn regarding the reaction mechanism and the nuclearity of the active catalyst.

Pd-Catalyzed Tsuji−Trost Allylic Substitution
Another example concerns the Pd-catalyzed Tsuji-Trost allylic substitutions and its combination with the V-catalyzed Meyer-Schuster rearrangement (Scheme 3) [21,22]. This coupling between propargylic alcohols and allylic carbonates promoted cooperatively by palladium and vanadium Tsuji-Trost O-allylations are generally considered to be catalyzed by mononuclear palladium complexes. Indeed, the reaction of Scheme 3 was initially modeled using a mononuclear complex in which the diphosphine ligand binds bidentately to the metal [23]. Although these calculations could establish several important aspects of the reaction mechanism, the results were not fully compatible with the experimental outcome. In particular, the dependence of the product distribution on the palladium loading could not be reproduced. On closer examination, it turned out that the particular diphosphine ligand used in the experiments, 1,1-bis(diphenylphosphino)methane (dppm), promotes the formation of dinuclear palladium complexes, while other diphosphine ligands in the same family, like for example 1,1-bis(diphenylphosphino)ethane (dppe), yield mononuclear complexes. Calculations on various Pd(0) complexes with dppm revealed that dinuclear complexes are much more stable than the mononuclear counterparts (Fig. 3). The ring strain of the four-membered palladacycle in the mononuclear complex is alleviated when the dppm ligands are bridging two Pd ions, resulting in structures that are as much as 25 kcal/mol lower in free energy. Starting from the complexes with a (dppm) 2 Pd 2 core, a reaction mechanism was then calculated that better matched the experimental findings. In particular, the calculations showed that only one of the palladium ions is active in the reaction, while the second one is a spectator that can be in oxidation state 0 or II. This results in two different catalytic cycles that have to be taken into account and that influence the final product distribution. In total, the calculations showed that the overall palladium-vanadium combined reaction consists of nine interconnected catalytic cycles. Using the obtained DFT energies, numerical kinetics simulations were then employed to work out the final reaction outcome and gain deeper insights into the factors governing the reactivity and selectivity of the reaction.
This case shows that a seemingly small variation of the ligand structure (dppm vs. dppe) can have a large impact on the reaction mechanism, and, again, a careful computational speciation study is crucial in order to model the reaction accurately and reproduce the experiments.

Cu-Catalyzed Cycloadditions
The final examples discussed here concern the Cu-mediated activation of terminal alkynes to form Cu-acetylide that then participates in cycloaddition reactions, such as the Huisgen cycloaddition to azides (the famed Click chemistry reaction, Fig. 3 Lowest energy di-and mononuclear Pd(0) complexes obtained from the computational speciation, with dppm as a ligand. Relative free energies (kcal/mol) are given in kcal/mol Scheme 4 a Cu-catalyzed azide-alkyne Huisgen cycloaddition, and b Kinugasa reaction 12 13 [Cu] R 1 H +
Coordination of the alkyne to a copper ion lowers its pK a significantly, making it possible for a base to abstract the proton and form a Cu-acetylide that can react further [24]. However, the calculations showed that this deprotonation, using triethylamine as a model base, is quite endergonic, by ca 9 kcal/mol [25]. Although this value is not very high, it has to be added to the barriers of the following steps and might therefore result in unfeasible overall barriers.
In the case of the Kinugasa reaction, the calculations could demonstrate that the presence of the second copper ion provides a considerable additional lowering of the pK a of alkyne. From being endergonic by 9 kcal/mol, the deprotonation by Et 3 N was found to be exergonic by 3 kcal/mol when two copper ions were employed (using phenanthroline as a copper ligand).
Comparison of the energy profiles of the initial steps of the Kinugasa reaction employing one or two copper ions is shown in Fig. 4. Apart from the further acidification of the alkyne, the second copper ion provides additional stabilization to the subsequent intermediates and transition states of the formal cycloaddition, which was shown to take place in a stepwise fashion. The second copper center has thus profound implications on the energetics of the reaction, and modeling the reaction using a mononuclear copper species does not yield feasible barriers.
In the case of the Cu-catalyzed azide-alkyne cycloaddition (CuAAC), both experimental and theoretical studies have demonstrated that two copper ions are involved in the reaction [26][27][28][29][30][31]. Although the precise reaction mechanism may vary depending on the reaction conditions, a general mechanism based on the previous mechanistic studies is summarized in Scheme 5. The alkyne reacts with two Cu(I) complexes to give Cu-acetylide 16. After coordination of azide to one of the copper ions, C-N bond formation takes place to produce a six-membered metallacycle 18, which undergoes a ring contraction and metal dissociation to give triazole intermediate 19. Subsequent protonation yields final product 13. Thus, also in this case, the second Cu ion facilitates both the deprotonation and that following ring formation step.

Conclusions
In this short account, we discuss a number of examples drawn from our own research where-much to our initial surprise-the active species in a catalyzed organometallic transformation turns out to be bimetallic. As discussed in the introduction, this turned out to be a surprise because the pre-catalyst species was monomeric and textbook mechanisms could be written with only a single metal center. In some cases, the intervention of multiple metals could be proved by experimental kinetics studies, while in others, the proposed mechanism is supported by it being better able to account for the observed reactivity and selectivity of the system. Hence, while the usual caveats This is one example of a general feature of contemporary computational studies: they by and large only provide information on those mechanisms which the computational chemist has decided to study. In that sense, the computations play the role of filling in the details and checking the possibility of a mechanism which was suggested based on prior examples or chemical intuition, rather than being 'predictive' in the strong sense of the term. As computational studies mature, methods become more accurate, and aspects such as computational speciation start to become systematically included, it is likely that the exploration of chemical reaction space needed to assess possible reaction mechanisms will allow an even greater degree of predictivity in studying catalysis. In the meantime, we hope that researchers in the field will remember that quite often, unity of catalyst metals makes strength-that two is better than one.