Ramiro F. Quijano-Quiñones*a,
Jareth Guadarrama-Morenoa,
Mariana Quesadas-Rojasbc,
Gonzalo J. Mena-Rejón
d,
Carolina S. Castro-Seguraa and
David Cáceres-Castillo
d
aLaboratorio de Química Teórica, Facultad de Química, Universidad Autónoma de Yucatán, Mérida, Yucatán, Mexico. E-mail: ramiro.quijano@correo.uady.mx
bPosgrado en Ciencias del Mar y Limnología, UNAM, Mexico
cEscuela Nacional de Educación Superior, UNAM, Mérida, Mexico
dLaboratorio de Química Farmaceútica, Facultad de Química, Universidad Autónoma de Yucatán, Mérida, Yucatán, Mexico
First published on 16th February 2021
Acrolein dimerization is a intriguing case since the reaction does not occur to form the electronically preferred regioisomeric adduct. Various explanations have been suggested to rationalize this experimental regioselectivity, however, none of these arguments had been convincing enough. In this work, the hetero Diels–Alder acrolein dimerization was theoretically investigated using DFT and MP2 methods. The influence of nucleophilic/electrophilic interactions and non-covalent interactions (NCI) in the regiospecificity of the reaction were analyzed. Our results show that the NCI at the transition state are the key factor controlling the regiospecificity in this reaction. Besides, we found that the choice of calculation method can have an effect on the prediction of the mechanism in the reaction, as all DFT methods forecast a one-step hetero Diels–Alder acrolein dimerization, while MP2 predicts a stepwise description for the lower energy reaction channel.
Traditionally, the endo rule has been explained mainly by SOIS, but, in recent years these interactions have been subject to controversy and NCI have been invoked to rationalize endo/exo stereoselectivities.6 Nevertheless, the influence of NCI on the regiochemical outcomes in DA and HDA, on the contrary to the one exerted by stereochemical factors, remains unclear. Since the first documented example of regiospecificity controlled by NCI in DA reactions between substituted acetylenes and tetrazines by Cioslowski et al. in 1993,3 few examples of DA and HDA regioselectivity mainly controlled by NCI have been reported. And many questions regarding the influence of NCI in regiocontrolling this type of reactions remain unanswered.
As is well known, acrolein thermally dimerizes through a regiospecific HDA cycloaddition (Scheme 1), producing exclusively isomer 27–10 and it is a textbook example of the regioselectivity of Diels–Alder reactions.11 The activation barrier energy for this reaction has been measured in a non-polar solvent (19.6 kcal mol−1, heptane) and is almost unaffected by the polarity of the media, indicating a non-polar mechanism. The proximity of the electron-withdrawing group to the unsaturated bond confers to acrolein strong electrophilicity,10,11 and the resonance forms of 1, indicate that the beta carbon is the most electrophilic atom, whereas the oxygen atom is the most nucleophilic. This is in complete agreement with the calculated reactivity index for this moiety.12 Therefore, according to electrophile–nucleophile interactions, adduct 3 should be the preferred regioisomeric product.
Until now, various explanations such as HOMO/LUMO interacting orbitals and SOIS has been adduced by different authors to rationalize the dimerization regioselectivity in acrolein dimerization.10,13–18 Even the regiospecific acrolein dimerization has been used as an example of the frontier orbital approach validity to regiochemical issues. But, due to the small differences in the coefficients of the HOMO/LUMO interacting orbitals,10,11 the FMO model could not explain the unique formation of the adduct 2. Thus, the origin of the regioselectivity in this reaction remains unclear.
In this work, we report the influence of the NCI in the regiospecificity of acrolein dimerization. To this end, we performed an accurate theoretical study on the NCI in the transition states imply in acrolein dimerization. The NCI was studied through the reduced gradient isosurfaces. The strength of the NCI was quantified by integrating their ρ(r) in defined ranges. The results show that the regiospecificity in the acrolein dimerization is being primarily controlled by the NCI. It is important to note that the analysis of the changes in electron density of the NCI along the reaction channels lies within the recent proposed Molecular Electron Density Theory (MEDT).19
Since the experimental activation energy9 was measured in an inert solvent and the experimental results show that the activation energy was almost unaffected by the solvent polarity, all calculation was made in the gas phase. Zero-point vibrational energy corrections were applied without scaling for all the examined structures.
The nature of all the stationary points was determined by diagonalizing the Hessian matrix at the same level. The reactants and products were identified from the vibrational analysis with all normal modes had real frequencies, and the TS with only one normal mode with imaginary frequency, corresponding to the movement in the direction of the reaction coordinate. Also, the intrinsic reaction coordinate (IRC) calculations were performed.42
The global electrophilicity index (ω) was calculated from the chemical potential (μ) and the chemical hardness (η) using the expression ω = (μ2/2η) according to Parr et al.43 The ω indicates the stabilization in energy when the chemical system acquires an additional charge from the surroundings. The η and the μ quantities can be approached in terms of the electron energy of the frontier molecular orbitals HOMO (EH) and LUMO (EL) as:
![]() | (1) |
η ≈ EL − EH | (2) |
In contrast to ω, several approaches have been established for the calculations of the nucleophilicity index (N). In this work, we have used the empirical approach to nucleophilicity proposed by Domingo et al.44 since there is evidence of its capability to predict the nucleophilic behaviour of organic molecules,45 in this scheme N is defined as
N = EHOMO(Nucleophile) − EHOMO(TCE) | (3) |
To characterize the most nucleophilic/electrophilic centers, we calculated the spin density on the radical anion of the electrophile and the radical cation of the nucleophile, using the Gaussian 09 (ref. 20) code, whereas the Spartan16 (ref. 46) software was used to generate the spin density maps. Then, the Parr function local reactivity indexes were calculated using the following equations.12
P+(r) = ρrcS(r) | (4) |
P−(r) = ρraS(r) | (5) |
Nk = NP−K |
ωk = ωP+K |
To study the influence of NCI in the activation energy of acrolein dimerization, calculations of the reduced density gradient was performed. In this method, the maximum variations in the contributions to the Laplacian, along with the axes, correspond to the eigenvalues (λi) of the electron-density Hessian matrix. The sign of λ2 enables us to distinguish between the different types of weak interaction (repulsive and attractive), while the electron density lets us assess the interaction strength. We depict the low gradient isosurfaces, over the range of −0.5 < Sig(λ2) ρ < 0.5. Such NCI calculations were performed using the NCIPLOT program,47 and the resulting isosurfaces were visualized with Visual Molecular Dynamics (VMD) software.48
Finally, he NBO analysis to calculate the global electron density transfer (GEDT) was performed with the version included in Gaussian 09.49–53
The eight possible TS are depicted in Fig. 1 where the label of each transition state is formed by combined all the symbols for the different combinations. To analyse the evolution of the potential energy surface for the dimerization in all reaction channels, we localize the eight TS (Fig. 1). Since results with all the functionals were qualitatively the same, only those obtained the results with M06-2X are presented, for simplicity.
![]() | ||
Fig. 1 Transition state structures for acrolein cycloaddition dimerization at the M06-2X level of theory. |
We calculated the distance of selected bonds (Scheme 2), the asynchronicity degrees (Δd), and the GEDT for all TS. The GEDT was computed by sharing the natural charges obtained from the NBO analysis.20,49
![]() | ||
Scheme 2 Nomenclature for selected bond lengths in the transition state structures for acrolein dimerization. |
The results, collected in Table 1, display that the beta carbon bond (d2) is the more developed in all the TS and show that all TS belongs to asynchronous process with the lowest energy TS 2CN presenting the higher asynchronicity degree (Δd). The GEDT at the most favourable TS 2CN (0.06e), allowed us to establish that acrolein dimerization as a non-polar reaction within the scale introduced by Domingo et al.54 Besides, has been established that with the increase of the polar character of DA reaction, the nucleophilic–electrophilic interactions at the TS are the mainly factor controlling reaction, whereas in a non-polar DA reaction this interactions cannot predict the regioselectivity.54 It is important to note, that the non-polarity of this reaction is in agreement with the experimental results. Finally, the global reactivity indexes of the acrolein, which show that this moiety is a good electrophile (ω = 1.6 eV) and good nucleophile (N = 1.5 eV).
d1 | d2 | Δd | GEDT | |
---|---|---|---|---|
2CN | 2.30 | 1.87 | 0.43 | 0.06 |
2TN | 2.17 | 2.00 | 0.17 | 0.09 |
2CX | 2.21 | 1.93 | 0.28 | 0.10 |
2TX | 2.16 | 1.98 | 0.18 | 0.10 |
3CN | 2.24 | 1.92 | 0.32 | 0.04 |
3TN | 2.13 | 2.01 | 0.12 | 0.07 |
3CX | 2.18 | 1.99 | 0.19 | 0.06 |
3TX | 2.14 | 2.01 | 0.13 | 0.06 |
In order to check the conclusions based on the GEDT values, we have calculated Parr functions (Fig. 2) which are in consensus with those calculated in the ref. 12 using a B3LYP/6-31G* method. The Parr functions indicate the most nucleophilic centre in the dienophile, is the alpha conjugated carbon atom (P−K = 0.15) whereas the beta conjugated carbon (P+K = 0.55) is the most electrophilic center in the diene. Therefore, as we expected from above results, the Parr indices cannot be used to predict the regioselectivity in this polar reaction. Thus, it is necessary to explore alternative factors that can explain the experimental observations, by analysing the transition states.
TS | MP2 | B3LYP | M06-2X | ωB97X-D | ωB97X | |||||
---|---|---|---|---|---|---|---|---|---|---|
Ea | Ereacc | Ea | Ereacc | Ea | Ereacc | Ea | Ereacc | Ea | Ereacc | |
2CN | 9.8* | −24.9 | 20.5 | −16.8 | 16.7 | −30.0 | 17.8 | −28.1 | 21.1 | −29.4 |
2CX | 14.5 | −25.2 | 22.3 | −16.2 | 19.1 | −30.4 | 20.2 | −28.2 | 23.3 | −29.4 |
2TN | 19.2 | −22.9 | 26.1 | −14.7 | 23.0 | −28.2 | 23.9 | −26.3 | 26.9 | −27.8 |
2TX | 19.9 | −23.2 | 26.8 | −14.1 | 23.8 | −28.6 | 24.6 | −26.4 | 27.6 | −27.8 |
3CN | 18.3 | −24.3 | 23.9 | −14.9 | 18.9 | −29.1 | 20.9 | −26.7 | 23.7 | −27.9 |
3CX | 20.9 | −23.1 | 25.0 | −14.9 | 20.9 | −28.0 | 22.6 | −26.3 | 25.4 | −27.1 |
3TN | 21.7 | −22.2 | 26.7 | −12.7 | 22.2 | −27.3 | 23.5 | −24.9 | 26.1 | −26.3 |
3TX | 22.0 | −21.1 | 26.6 | −12.8 | 22.6 | −26.1 | 24.0 | −24.4 | 26.6 | −25.5 |
It is known that the B3LYP functional overestimated the Ea in DA and HDA reactions where the number of double bonds is reduced, since this functional cannot provide a balanced description of reaction that involve double to single bond transitions, as it energetically favours π-orbitals over σ-orbitals.26,66,67 For these reasons, in the acrolein dimerization case, B3LYP predicts reaction energies less exothermic than all the other methods, as we can check in Table 2. However, surprisingly, the B3LYP functional provide the best fit between the calculated activation energy and the experimental activation energy. Considering the factors discussed above, this B3LYP's good performance could be due to the well-known error cancellation in this functional68,69 rather than a good description of the electronic density on this system.
It is also notable that MP2 does not predict a one-step HDA reaction for the lower activation energy 2CN channel, as all attempts to find one-step pathways were unsuccessful. Instead, it predicts a stepwise reaction where the rate-limiting step of the mechanism is the first step with a Ea = 9.8 kcal mol−1, in clear disagreement with the experimental activation energy.8,9 In contrast, each one of the TS calculated with DFT methods belong to a one-step HDA cycloaddition reaction. Moreover, all DFT methods show that cisoid conformation lead to more stable products and transition states, compared to transoid conformations, whereas this result is not predicted by MP2, neither in the activation energies nor in the transition states. The MP2's failure in the description of this reaction suggests that any MP2 prediction about DA and HDA mechanism should be taken carefully, and we suggest the use of another method to confirm the results for the MP2's geometry optimization of the TS.
The ωb97X functional overestimated the activation energy of the all TS respect to the ωb97X-D results by 2.3 (3TX) up to 3.3 kcal mol−1 (2CN). As the ωb97X functional is a version of the ωb97X-D functional without empirical dispersion correction, the differences between the results of both functionals are likely caused by London dispersion interactions. These results indicate that the NCI may play an important role in the stabilization of the transition states for acrolein dimerization, and could be the explanation of the regiospecificity of this reaction.
In this sense, the quantification of the NCI is of crucial importance for understanding their importance in the chemical structure and reactivity. The integration of ρ(r) is defined range allows an estimation of the strength of the NCI70,71 and recently, the relationship between the NCI integrals and interaction between two fragments has been addressed.72,73 Therefore to analyse how NCI are involved in the stabilization of the transition states, a set of calculations have been performed integrating the electronic density of the NCI in the range corresponding to attractive (−0.1 < sign(λ2)ρ(r) < −0.02) and vdW interactions (−0.02 < sign(λ2)ρ(r) < 0.02) using the code NCIPLOT4.47 As results with all the functionals were qualitatively the same, for simplicity only the results with M06-2X are presented.
Fig. 3, shows the relation between the sum of the vdW and attractive integrals and the activation energy of all transition states and show that the NCI correctly account for the preference of the cisoid over the transoid conformations. Moreover, the energetic preference of all the TS is explained adequately by the NCI, as there is a significant linear correlation (R = −0.94, p = 0.003) between the NCI and the activation barriers.
Accordingly to the above outcomes, a NCI analysis through the reduced gradient density was performed on the transition state structures calculated at the M06-2X method, to establish the favourable NCI appearing between both fragments. By analysing the low gradient isosurfaces (Fig. 4), we can observe cisoid/endo transition states (2CN and 3CN) have more covalent interactions, respect to exo and transoid transition states. This is likely due to their more compact shape that allows more direct interaction between both fragments.
Furthermore, the NCI are more pronounced in TS 2CN, as the oxygens of both carbonyl groups can interact with the carbon of the carbonyl from the other acrolein moiety, forming a carbonyl–carbonyl interaction (CO⋯CO).74 This interaction accounts for the stabilization of 2CN respect to other reaction channels.
Even though in exo approximation the atoms are further away than in endo approximation, in the cisoid/exo transition states (2CX and 3CX) a weak hydrogen bond (HB) is formed with a bond length and bond angle of 2.56 Å and 112.8° for 3CX; and 2.49 Å and 131.1° for 2CX. This HB explains the energy stabilization of these transition states, respect to the transoid TS, where neither hydrogen bonding, nor the oxygen (dienophile)–carbon (diene) interactions are possible.
In summary, the analysis of the changes in electron density of the NCI along the reaction channels show that the regiospecificity observed in acrolein dimerization can be explained by NCI analysis of the TS, and clarifies the theoretical endo/exo preferences. Finally, all the NCI calculated in this work belong to regions with low density and low reduced gradient but with no zero gradients, therefore their study was not possible using an atom-in-molecules (AIM) approach.75,76
Footnote |
† Electronic supplementary information (ESI) available. See DOI: 10.1039/d0ra10084f |
This journal is © The Royal Society of Chemistry 2021 |