Sarote
Boonseng
,
Gavin W.
Roffe
,
John
Spencer
and
Hazel
Cox
*
Department of Chemistry, School of Life Sciences, University of Sussex, Falmer, Brighton, East Sussex, BN1 9QJ, UK. E-mail: h.cox@sussex.ac.uk
First published on 18th March 2015
The accuracy of DFT-optimised geometries of the symmetrical pincer palladacycles PdNCN and PdSCS, [ClPd{2,6-(Me2NCH2)2C6H3}] and [ClPd{2,6-(MeSCH2)2C6H3}] respectively, has been evaluated by investigating the performance of eight commonly used density functionals with four combinations of basis set, in reproducing their X-ray crystal structures. It was found that whilst the ωB97XD functional performed best over all, the PBE and TPSS functionals performed best when considering the palladium coordination geometry. The role of the donor atom in the stability and reactivity of the symmetric palladacycles, PdYCY, Y = N, S, or P, has been determined using Bader's Atoms in Molecules method to elucidate the nature of the bonding, and using a model formation reaction, which involves the C–H activation of the pincer ligand YCY by PdCl2. The calculations reveal distinct differences in the bond strength and nature of the interaction of Pd with the donor atoms Y, which support differences in the thermodynamic stability of the palladacycles.
The structure and bonding of palladacycles can be investigated using density functional theory (DFT), however an appropriate choice of functional, basis set, and relativistic effective core potential (ECP) for Pd, is essential for obtaining accurate structures, and hence thermodynamic and kinetic data.14 Several papers provide useful benchmarking studies of DFT geometries compared to experimental data including the recent study by Minenkov et al.15 The authors compared DFT-optimised geometries with experimental crystal structures for 18 ruthenium complexes, and a further 10 transition metal complexes. The results showed that the use of functionals that account for dispersion decrease the statistical error between experiment and theory compared to those that do not include dispersion corrections, with the ωB97XD functional providing the best overall results. However for accuracy around the metal centre the PBE and TPSS functionals also performed very well, and for organic ligands B3LYP performed very well. Waller et al.16 investigated the ability of 15 functionals to describe the optimised geometries of 19 second-row transition metal complexes compared to experimental gas-phase data, showing that hybrid functionals such as the PBE hybrid, B3PW91 and B3P86 provide the most accurate geometries. Clearly the choice of functional is dependent on the chemical nature of the structures studied.
The purpose of this paper is two-fold. First, to determine an optimum DFT methodology for the study of pincer palladacycles that combines accuracy with computational speed. This will be achieved by evaluating the performance of a range of functionals in their ability to reproduce the structural features of two experimentally characterised symmetric pincer palladacycles (Fig. 1). Second, to determine the role of the donor atoms, Y = N, S or P, in the symmetric palladacycles, PdYCY. To compare stability and reactivity, the steps involved in a simple formation reaction, which involves the C–H activation of the pincer ligand YCY by PdCl2, will be calculated, and the strength and nature of the Pd–L bonding in PdYCY, where L = Y, C, Cl, evaluated.
The optimum methodology for geometry optimisation in the present work was found to be PBE/6-31+G(d,p)[SDD] and was therefore used to calculate the stationary points along the formation reaction pathway of PdNCN and PdSCS (Fig. 1) along with the phosphorus analogue PdPCP ([ClPd{2,6-(Me2PCH2)2C6H3}]) (Fig. 2).
All minima and transition states were verified by, respectively, the absence or presence of a single imaginary mode. Single point energy calculations were performed at the ωB97XD/6-311++G(2df,2p)[SDD] level of theory, as the ωB97XD functional25 includes important non-covalent interactions and has been shown to provide accurate energetics.32–35 The zero-point energy or the Gibbs free energy corrections obtained at PBE/6-31+G(d,p)[SDD] level of theory were applied to the single point energy calculations unscaled.
A topological analysis of the electron density was performed using the Atoms in Molecules (AIM) method,36,37 as implemented in the Multiwfn program,38 to analyse the bonding using the ωB97XD/6-311++G(2df,2p)[DGDZVP]//PBE/6-31+G(d,p)[SDD] model chemistry. The use of ECPs are not recommended for AIM analysis as bond paths cannot be traced,39 hence the relativistic DGDZVP all electron basis set was used to treat palladium.40
The results from the Quatfit analysis when all atoms have equal weighting show that there is very little difference between the performances of the functionals (Fig. 3), with the RMS values ranging from 0.067 to 0.080 Å.
Fig. 3 RMS error for optimised PdNCN and PdSCS structures compared to X-ray crystal structures using Quatfit program with equal weighting for all non-hydrogen atoms. |
The Quatfit analysis using only palladium and directly bound atoms (Fig. 4), shows that these bonds are better reproduced and show a slightly greater variation in the errors, with the GGAs PBE and BP86, along with the meta-GGA TPSS performing best. However, the difference between the best and worst performing functional is just 0.019 Å.
Fig. 4 RMS error for optimised PdNCN and PdSCS structures compared to X-ray crystal structures using Quatfit program for Pd–L bonds. |
When the basis set choice is studied, it is very clear that the SDD ECP outperforms LanL2DZ (Fig. 3 and 4). There is negligible difference between the 6-31G(d) and 6-31+G(d,p) basis sets, however the inclusion of the diffuse functions is desirable because this provides greater flexibility without a substantial increase in computational time. Therefore, the 6-31+G(d,p)[SDD] basis set was chosen as the optimum basis set for geometry optimisation.
Minenkov et al. showed that functionals that do not account for dispersion interactions systematically overestimated internuclear distances, whereas when dispersion was included these overestimations were matched with underestimates, resulting in very small mean signed errors (MSE).15 In the study by Waller et al. it was also found that standard DFT functionals overestimated bond distances.16 Therefore complementary to the Cartesian coordinate approach utilised by the Quatfit analysis, in order to gain an insight into any systematic over or underestimation of bond distances and interatomic distances, a total of 422 interatomic distances and 34 chemically bonded distances calculated using the 6-31+G(d,p)[SDD] basis set and the eight density functionals were compared with the experimentally determined distances. The mean unsigned error (MUE): the average of the absolute deviation between calculation and experiment, and the mean signed error (MSE): the average deviation, were calculated (Fig. 5). The MSE results indicate functionals containing dispersion: M06, M06L and ωB97XD, with the exception of B97D, experience a greater degree of fortuitous cancellation of errors by the cancellation of some of the overestimated distances with underestimated bond distances as reported by Minenkov et al.,15 for both the bonded and interatomic distances. There are also slightly lower MUE values for functionals accounting for dispersion (except B97D), however once again the differences between all the functionals were relatively small with MUE values for bonded distances between 0.011 and 0.022 Å, and for all interatomic distances slightly larger between 0.022 and 0.046 Å.
Fig. 5 Average bond and interatomic distance MSE and MUE (excluding hydrogens) for PdNCN and PdSCS compared to their X-ray crystal structures using 6-31+G(d,p)[SDD] basis set. |
Of key importance is the Pd–L bonding environment. Fig. 6 depicts the MUE and MSE values for the eight Pd–L distances of the PdNCN and PdSCS structures. It is clear that, as per the Quatfit data, the PBE and TPSS functionals reproduce the metal ligand environment extremely well.
Fig. 6 Comparison MUE and MSE of Pd–L bonds of both PdNCN and PdSCS using the 6-31+G(d,p)[SDD] basis set. |
Taking into account the data from the MUE and MSE of all bond and interatomic distances, all the DFT functionals considered predict expanded structures, but the errors are reasonably small with ωB97XD having the smallest error. The non-hybrid functionals PBE, BP86 and TPSS provide the best results for the Pd–L environment, and provide advantageous calculation time compared to the other functionals which appear higher on Jacob's ladder.41 The importance of the Pd–L environment, and the desire for an accurate but computationally efficient method for optimisation of large molecular structures led to the choice of the PBE functional, although the use of TPSS or ωB97XD would be equally appropriate. To confirm this, the Gibbs free energy of formation of the palladacycles via PdCl2 using both the PBE- and ωB97XD-optimised geometries were compared and found to differ by less than 2 kJ mol−1 (see Table 1). Therefore the methodology of choice for optimisation in the proceeding section is PBE/6-31+G(d,p)[SDD]. Single point energies using ωB97XD/6-311++G(2df,2p)[SDD] at the optimised geometry were performed to ensure accurate energetics as it has been shown that the inclusion of non-covalent interactions is essential for accurate thermodynamics and kinetics, for example, it has been shown that they constitute a significant proportion of the binding energy in various transition metal phosphine complexes.32,33
Complex | ΔG0/kJ mol−1 | |
---|---|---|
PBE | ωB97XD | |
PdNCN | −207.3 | −206.6 |
PdSCS | −213.2 | −214.9 |
PdPCP | −318.7 | −320.4 |
The data in Table 1 show that the PdPCP complex is the most thermodynamically stable, and that the formation of all three complexes is spontaneous (ΔG0 < 0).
The simplest concerted formation reaction,42,43 was studied for PdNCN, PdSCS and PdPCP (Scheme 2) with the energy profile for each reaction provided in Fig. 7. The purpose of this study was to gain insight into the reactivity of the pincer ligands and to provide a reference point when comparing symmetrical and unsymmetrical pincer palladacycle structures. This mechanism involves several steps. Initially the metal coordinates to one of the donor atoms of the pincer ligand and then inserts via a concerted C–H bond activation to form a coordinately-unsaturated Pd centre weakly bound to HCl which is eliminated before the second donor atom of the ligand coordinates to the metal centre to form the PdYCY palladacycle.
Scheme 2 An unassisted formation reaction pathway for PdNCN, PdSCS and PdPCP using monomeric PdCl2. Int = Intermediate and TS = Transition State. |
Fig. 7 ZPE corrected energy profiles for the formation reaction pathway (Scheme 2) of PdNCN, PdSCS and PdPCP. |
The concerted C–H bond activation step proceeds viaTS1. C–H activation is one of the most widely studied and important topics in organometallic chemistry due to its importance in catalysis.44,45 The C–H bond activation energies for PdNCN, PdSCS, and PdPCP are 118.0, 101.5 and 92.5 kJ mol−1, respectively (Fig. 7). The observation that the barrier to C–H activation in this concerted activation pathway is lowest for PdPCP is in agreement with a similar investigation into the C–H activation of CH4 by late transition metal pincer complexes.45 This step has the highest activation energy barrier in the mechanism and thus constitutes the rate-determining step. The second activation step consists of the coordination of the second donor atom to the palladium centre viaTS2. This occurs via the inversion of the second donor atom into the position to begin coordination to the metal centre for PdNCN and PdPCP, and the rotation of the ligand arm for PdSCS. Due to the different mechanisms, the TS2 barrier for PdNCN and PdPCP is considerably larger (34.9 and 87.4 kJ mol−1, respectively) than that of PdSCS (5.4 kJ mol−1). Decoordination of one of the donor atoms from the Pd centre is likely to be a key step in the catalyst activation.46
The AIM parameters at the key BCPs are provided in Table 2. The magnitude of the electron density ρ(r) at the BCP indicates the strength of the interaction. The Laplacian of the electron density ∇2ρ(r) indicates the nature of the interaction; ∇2ρ(r) < 0 indicates regions of local electronic charge concentration corresponding to a covalent-type interaction, whereas ∇2ρ(r) > 0 indicates density depletion corresponding to an ionic or closed-shell interaction.39 Additionally, the sum of the gradient kinetic electron density G(r) and the potential energy density V(r), provides the total electron energy density H(r), which further elucidates the nature of the bond.
Structure | Bond type | ρ(r) | ∇2ρ(r) | H(r) | ε | δ(Pd–L) |
---|---|---|---|---|---|---|
PdNCN | Pd–N1 | 0.087 | 0.383 | −0.014 | 0.033 | 0.908 |
Pd–N2 | 0.086 | 0.383 | −0.014 | 0.033 | 0.908 | |
Pd–C | 0.156 | 0.236 | −0.070 | 0.081 | 1.126 | |
Pd–Cl | 0.064 | 0.238 | −0.011 | 0.198 | 1.219 | |
PdSCS | Pd–S1 | 0.091 | 0.257 | −0.025 | 0.040 | 1.157 |
Pd–S2 | 0.091 | 0.257 | −0.025 | 0.040 | 1.157 | |
Pd–C | 0.137 | 0.223 | −0.055 | 0.061 | 1.051 | |
Pd–Cl | 0.068 | 0.237 | −0.013 | 0.114 | 1.234 | |
PdPCP | Pd–P1 | 0.101 | 0.183 | −0.036 | 0.022 | 1.124 |
Pd–P2 | 0.101 | 0.183 | −0.036 | 0.022 | 1.124 | |
Pd–C | 0.126 | 0.224 | −0.047 | 0.050 | 1.029 | |
Pd–Cl | 0.069 | 0.232 | −0.013 | 0.085 | 1.258 |
The data in Table 2 indicate that, for all PdYCY complexes the Pd–C bond is the strongest, indicated by the largest ρ(r), and the Pd–Cl bond is the weakest. Furthermore, the Pd–P bonds are stronger than the analogous Pd–S and Pd–N bonds at the BCP. This is corroborated by the fact that PdPCP is the most thermodynamically stable (Table 1), and the comparatively larger energetic gains of coordination of P to Pd in the formation of Int1 and Product (Fig. 7). The bond strengths of the Pd–N and Pd–S bonds in PdNCN and PdSCS respectively are similar, with slightly stronger Pd–S interactions.
In Suzuki–Miyaura catalysis,47 and other types of catalysis where the active species are Pd(0) species, the active species must be generated from the palladacycle precatalysts. A necessary step in this activation is the reductive elimination of Pd from the Pd–C bond.46 In an investigation into the reductive elimination of ethane from Pd(PR3)2(Me)2 by Sajith and Suresh it was found that the nature of the bond to be eliminated is more important than their strengths.40 They found that a more ionic bond results in a higher activation barrier for the reductive elimination whereas a more covalent bond is easier to cleave. Therefore in the present work the nature of the bonds present in PdNCN, PdSCS and PdPCP, is determined, as this may prove important in fine-tuning the reactivity of pincer palladacycles.
Using the combination of ∇2ρ(r), H(r)40,48 and the delocalisation index δ(Pd–L) which determines the degree of covalency,49,50 the nature of the bonding between palladium and the donor atoms can be determined. All the Pd–L bonds exhibit partial ionic and covalent character, the former indicated by the positive ∇2ρ(r) and the latter by the negative H(r). This type of character has been described as “transit closed shell” bonding, with some covalent contribution.51,52 These transit closed shell interactions are typical in many transition metal donor atom bonds.40,50,53–55 The extent of sharing of electrons between the bonds, δ(Pd–L), provides an independent measure of the degree of covalency. The δ(Pd–L) values range from 0.9–1.3 indicating that all of the Pd–L bonds have significant covalent character.
The Pd–N bonds have most ionic character, demonstrated by the largest ∇2ρ(r) and H(r), and smallest δ(Pd–L). The values of ∇2ρ(r) and H(r) at the BCPs indicate that the Pd–P bonds have the most covalent character, however the Pd–S bonds have slightly larger delocalisation indices (Table 2). The bond degree parameter |H(r)/ρ(r)|, which can be interpreted as the total pressure per electron density unit,51 also indicates that the Pd–P bond is the strongest (|H(r)/ρ(r)| = 0.356) and the Pd–N the weakest (|H(r)/ρ(r)| = 0.161).
The bond ellipticity, ε, can be used to determine the degree of σ and π character in the Pd–L bonds. It measures the extent to which the electron density is unequally distorted away from the bond axis.56 The low ε of all the Pd–L bonds reflect their predominant σ character, with the Pd–Cl interactions exhibiting a greater π character contribution.
It is clear from this analysis that all of the Pd–L bonds are characterised by low ρ(r), positive ∇2ρ(r) and negative H(r) indicating partial ionic and covalent character. However, it is also clear that there are subtle differences between the Pd–L interactions; the relative magnitude of the parameters indicates that the Pd–P bond is most covalent and the Pd–N bond most ionic. This varying nature of the bonds in each palladacycle could potentially have implications in reactivity.
The role of the donor atoms in the stability and reactivity of the symmetric pincer palladacycles, PdYCY, Y = N, S or P, was evaluated. The simple formation reaction of the YCY pincer ligand with PdCl2, via C–H activation, was determined in order to give a measure of the relative stability and reactivity as a function of the donor atom in the absence of other effects (such as intra- or inter-molecular assisted formation and solvent effects). It was found that both thermodynamically and kinetically PdPCP was the most stable and PdNCN was the least stable. This was supported by a topological analysis of the electron density in the PdYCY complexes. The data suggested that Pd–P bonds were the strongest and the Pd–N bonds the weakest. Furthermore, it was found that the Pd–P bond had more covalent character than the Pd–S and that the Pd–N bond had the most ionic character.
The calculations clearly revealed distinct differences in the bond strength and nature of the interaction of Pd with the donor atoms N, S and P, which supported the differences in the stability, namely that the PdPCP was the most thermodynamically stable (in the absence of external effects). The varying nature of the bonds in each palladacycle could potentially have implications in their reactivity. Exploitation of these differences by considering unsymmetric YCY’ complexes, both experimentally and theoretically, are underway.
Footnote |
† Electronic supplementary information (ESI) available: Energies and Cartesian coordinates of the stationary points along the formation reaction pathway are provided. See DOI: 10.1039/c5dt00031a |
This journal is © The Royal Society of Chemistry 2015 |