Exploring CO₂ activation mechanisms with triphenylphosphine derivatives: insights from energy decomposition and deformation density analyses

Abstract

This study focuses on the reaction mechanisms involving triphenylphosphine (PPh₃) derivatives, benzyne, and CO₂, giving mechanistic insights into two competing pathways: Path a, which involves direct C–P bond formation, and Path b, which progresses via a [2 + 2] cycloaddition. Comprehensive computational analysis by energy decomposition analysis (EDA) and deformation density insights was employed to elucidate the electronic and steric factors influencing the reactivity and selectivity of PPh₃ derivatives. The results reveal that Path b is energetically and kinetically favored. In Path a, substantial repulsive interactions (ΔE_rep), especially for electron-withdrawing substituents, hinder C–P bond formation, making this pathway unfavorable, while Path b benefits from compensatory effects between interaction energies, with electron-releasing para-substituents, such as NHMe and OMe, increasing stabilization by enhancing ΔE_orb contributions. Substituents in meta positions show greater distortion energies (ΔE_dist), which limit their stabilizing effects compared to para-substituents. The deformation density analysis of transition states (TS1(b) and TS2(b)) emphasizes the crucial role of Pauli deformation (Δρ^Pauli) and orbital deformation (Δρ^Orb) in modulating stability. Para-substituents exhibit stronger electronic effects, reducing ΔE_int more effectively than meta-substituents, which increase ΔE_dist. This positional dependence underscores the importance of substituent design in optimizing reactivity.

---

1. Introduction

Efficiency in the development of reaction mechanisms with the participation of triphenylphosphine derivatives remains one of the main goals of organic chemistry owing to their versatility in promoting complex transformations. Energy Decomposition Analysis (EDA) has been instrumental in elucidating key interactions in transition states and intermediates, providing a deeper understanding of selectivity in organic reactions.^1–5 EDA is widely used in organic chemistry to break down the total interaction energy between reacting species into several components, such as electrostatic, exchange, polarization, and charge transfer contributions.⁶ These methods quantize the individual contributions of factors affecting reactivity, reaction mechanisms, and product formation.⁷ Thus, the forces operating are more clearly understood. Most importantly, EDA allows the differentiation between various interaction types that influence both the stability of the reaction's intermediates and transition states.⁸ A possible example in that direction is catalysis: EDA can identify crucial interactions between catalysts and substrates useful in finding optimal reaction conditions or in the design of more efficient catalytic systems.^9–11 In any case, EDA holds a special position in the investigation of electron transfer processes, since it is able to dissect the role of charge transfer in the case of reactions involving radical or ionic intermediates. The EDA methods have also been applied to the analysis of other reaction types, including noncovalent interactions, hydrogen bonding, and π–π stacking in organic synthesis. The capability to visually quantify such interactions has become important in designing new organic materials, including catalysts and molecular devices. Because of this, EDA has become a very integral part of mechanistic studies in organic chemistry due to the view it presents about the molecular nature of interactions and its power as an approach to improvements in efficiency and selectivity in chemical reactions.^12–16

Understanding the mechanisms of CO₂ transformations is not only important for the development of green chemistry but also for organic chemists to develop efficient strategies for carbon incorporation. Mechanistic studies in CO₂ transformation, play a critical role in environmental sustainability, economic viability, and the advancement of green chemistry.^17,18 For further details, researchers often look into mechanistic research using tools like Density Functional Theory (DFT) calculations to model and predict these catalytic processes, driving progress in this field. Despite its broad utility, EDA has been underutilized in understanding the mechanistic pathways of CO₂ transformations, especially in the context of triphenylphosphine derivatives. This work tries to fill this gap by using EDA to provide mechanistic insights and guide substituent design. However, a literature review shows that in some limited investigations, the EDA has been applied very successfully. Notably, the absolutely localized molecular orbital (ALMO)-EDA approach has been applied to evaluate electronic interactions in catalytic systems involving CO₂, examining factors like polarization, charge transfer, and dispersion that contribute to the reaction's efficiency.¹² Also, another study utilizes the ALMO-EDA (solv) approach, incorporating solvation effects to analyze intermolecular interactions in CO₂ reduction catalysts. This method allows researchers to separate interaction energies into components—such as electrostatics, polarization, and charge transfer—thereby offering insights into catalyst behavior in solution and the stabilization of intermediates involved in CO₂ reduction.¹⁹

Beyond CO₂ transformation, triphenylphosphine derivatives are critical in several organic transformations like Wittig reactions and cross-coupling, thereby underlining their importance in the design of reactions and catalysis.²⁰ Their role in developing more sustainable chemistry processes is increasingly relevant as researchers look for efficient ways to convert CO₂ into commercially valuable compounds.^21–24

Lin He and coworkers reported the synthesis of some zwitterionic phosphonium salts via the reaction of benzynes, some phosphine derivatives, and CO₂.²⁴ They proposed two paths (a, b) for the reaction, and regarding the experimental facts, they reported that benzyne first undergoes a [2 + 2] reaction with CO₂ to generate an intermediate having a four-membered ring. Then owing to the high ring strain of the four-membered ring, the intermediate undergoes a ring-opening reaction by a nucleophilic attack of triphenylphosphines. The experimental observations in this report were an inspiration for us to study the reported mechanistic paths and also the substituent effects on the mechanism. Thus, in this study, we provide a comprehensive computational analysis of the reaction mechanisms involving PPh₃ derivatives, including electron-releasing and withdrawing group, benzyne, and CO₂, with a focus on two competing pathways: one involving direct C–P bond formation and the other progressing through a [2 + 2] cycloaddition intermediate. By leveraging EDA and deformation density methodologies, we systematically investigate the influence of substituent effects on transition state stability and intermediate formation. Notably, we explore how electron-withdrawing and electron-donating groups impact repulsive (ΔE_rep) and attractive (ΔE_orb, ΔE_els) interaction energies, offering insights into the electronic and steric factors that govern the reactivity and selectivity of PPh₃ derivatives in complex organic reactions. The results not only advance the mechanistic understanding of PPh₃-based transformations but also suggest strategic modifications to enhance reaction efficiency and selectivity. Furthermore, the results illustrate the ability of EDA and deformation density methodologies in studying CO₂ transformation mechanisms. This work provides a valuable foundation for future studies aimed at designing PPh₃ derivatives tailored for specific organic transformations, potentially broadening the scope of PPh₃ applications in CO₂ transformations.

2. Computational details

Several methods have been worked out for approximating the exchange-correlation functional in Kohn–Sham DFT, including a division into local and nonlocal classes. Local spin-density approximation means a functional that depends only on the local values of the spin densities. However, LSDA is not sufficiently accurate for many applications, and the inclusion of some fraction of orbital-dependent nonlocal Hartree–Fock exchange led to the so-called hybrid functionals, including M06 and M06-2X.²⁵ These hybrid functionals are distinguished by the % Hartree–Fock exchange, denoted X here. Several previous studies appeared to demonstrate that this family of hybrid functionals, M06, had better performance than the B3LYP in most applications involving main-group thermochemistry, energy barriers for reactions, and noncovalent interactions. The M06 family has been applied to both closed-shell and open-shell systems, accounting for the ground-state spin quantum number as well as excited states of atoms.^26,27 The geometries of reactants, TS, and products were optimized by M06-2X/def2SVP^27,28 using the Gaussian 09 package.²⁹ Frequency calculations were performed on all structures at 298 K to confirm the minimum of the potential energy surface. The TS structures were obtained using Schlegel's synchronous transit-guided quasi-Newton (STQN) method to ensure that just one imaginary frequency existed for the calculated TSs. Moreover, the Intrinsic Reaction Coordinate (IRC) has been done to confirm the TSs are correctly assigned.³⁰ The solvent effects were evaluated in tetrahydrofuran (THF) as the solvent by a conductor-like polarizable continuum model (CPCM).³¹ Single-point energy calculations have also been performed at the M06-2X/def2tzvp level of theory for further improvement in the accuracy of the results.

Generally, the activation energy (ΔE^≠) of the studied transition states (TSs) could be decomposed into the distortion energy (ΔE_dist) of the involved fragments and the interaction energy (ΔE_int) between them (eqn (1)) using the distortion–interaction model.^6,32

Based on the energy decomposition analysis (EDA) method using dispersion-corrected DFT strategy (or the so-called sobEDA),³³ the interaction energy (ΔE_int) between the involved segments was divided into six terms including electrostatic (ΔE_els), exchange-reciprocal (ΔE_x + ΔE_rep), orbital interaction (ΔE_orb) and the sum of the DFT correlation energy (ΔE_DTFc) and dispersion interaction (ΔE_dc) reflects coulomb correlation (ΔE_c) terms (eqn (2)). Thus, with respect to eqn (3), the overall activation energy (ΔE^≠) is composed of the portions of steric effects (ΔE_steric), electronic effects (ΔE_elec), and Coulomb correlation effects (ΔE_c) in the transition state.

ΔE_int = ΔE_c + ΔE_rep + ΔE_els + ΔE_orb + ΔE_DFTc + ΔE_dc (2)

The other concept utilized for an explanation of the observed energies in transition states is deformation density.^34,35 This concept describes the shifts in electron density that occur when atoms interact to form a molecule, capturing the redistribution of electron density that distinguishes bonded atoms from those that are separate. Deformation density (Δρ^Total) can be decomposed into various components, with Pauli deformation density (Δρ^Pauli) and orbital deformation density (Δρ^Orb) being among the most informative (eqn (4)).

Δρ^Total = Δρ^Orb + Δρ^Pauli (4)

Δρ^Pauli is due to electron repulsion in terms of kinetic energy pressure and steric exclusion, while Δρ^Orb causes orbital relaxation, mixing, and consequent electron density redistribution. Fakhraee and Azami³⁶ specified kinetic energy pressure as a measurable form of steric repulsion, while Frenking et al.³⁷ associated Δρ^Pauli with like-spin electron repulsive interactions and equated Δρ^Orb with orbital mixing which leads to bonding rearrangements. Together, the reported findings underpin a mechanistic picture wherein the interaction between Δρ^Pauli and Δρ^Orb determines repulsive and attractive forces governing chemical bonding and material properties. The two-faced perspective supports analyses ranging from interpreting metal stacking fault energies to describing intermolecular interactions in molecular clusters.

These components are highly valuable when analyzed through the Natural Orbitals for Chemical Valence (NOCV) theory.^38,39 NOCV theory systematically dissects electron density to reveal bonding interactions, constructing natural orbitals that optimally represent chemical bonding. In this framework, Δρ^Pauli arises primarily due to constraints from the Pauli exclusion principle, showing electron density redistribution as electrons with the same spin repel each other to avoid overlap, thereby preserving spatial separation. Meanwhile, Δρ^Orb accounts for the changes in electron density due to the interaction of orbitals as a result of bonding and electron delocalization. This shows those areas of space where electrons have either concentrated or depleted due to the formation of bonds or simple overlap, whether a sigma or pi bond it is. Deformation density components allow for a more instructive view on the nature of molecular bonding both qualitatively and quantitatively. The research within the density functional theory with the use of NOCV theory and deformation density decomposition has been developed in understanding specific interactions affecting the stability of a molecule and bonding energy accordingly.⁴⁰ EDA analysis and deformation density components were evaluated using the MultiWFN 3.8.⁴¹

3. Results and discussion

Fig. 1 shows the studied triphenylphosphine derivatives in the reaction. Fig. 2 depicts both considered mechanisms a and b in which solid and dashed lines depict the progress of the reaction for para and meta substitutions, respectively. Fig. 3 and 1(S)† depict the corresponding Potential Energy Diagrams (PED) for mechanism b and a, respectively, in which the total Gibbs energies of CO₂, PPh₃, and benzyne are considered as the zero state, thus the energies of the other species are abstracted from the zero state.

Fig. 1 The optimized structure of PPh₃ derivatives; and the applied abbreviation during the studies in bolded form in parentheses.

Fig. 2 The studied mechanisms a and b for the reaction.

Fig. 3 The PED of the mechanism b.

Path a progresses through forming a complex (complex(a)) of benzyne structure and triphenylphosphine derivatives, then, zwitterion In1(a) (depicted in solid rectangular boxes in the PED of Mechanism a), suffering a negatively charged carbon atom, is formed via forming a C–P bond in Step1(a). Step2(a) includes developing a C–C bond between the carbon atom of CO₂ and the C2 atom of In1(a). None of these steps in Path a do not pass transition states, because they go through an electrostatic interaction between charged or partially charged atoms. Thus, our attempts to find corresponding transition states were not successful. Path b initiates by forming a complex between benzyne and CO₂ (complex(b)). Step1(b) is a [2 + 2] cycloaddition reaction which leads to In1(b) (depicted in dotted rectangular boxes in the PED of Mechanism b) after descending of TS1(b) (depicted in a solid rectangular box in the PED of Mechanism b). In Step2(b) the apt four-membered ring of In1(b) passes through TS2(b) (depicted in dashed rectangular boxes in the PED of Mechanism b) by nucleophilic attack of the phosphorous atom of PPh₃ to C1 atom, the final product is the outcome of this step. Table 1(S)† shows the calculated thermodynamic and kinetic parameters for the mechanisms.

3.1 HOMO/LUMO band gaps

Based on the Frontier Orbital Theory (FMO), the reaction can be considered as an interaction between high occupied orbital (HOMO) of PR₃ derivatives and low unoccupieolecular orbital (LUMO) of CO₂, benzyne or In(1)b. Table 2(S) shows the calculated HOMO/LUMO values for the involved species (in a.u.). In this table the energy gaps between are useful criteria to find whether the reaction progresses through Path a or b. In fact, the absolute value of the energy gap is lower for all derivatives. On the other hand, the values of the energy gap of is higher than the two others. Thus, it is unlike the reaction going through a reaction between PR₃ derivatives and CO₂. Also, the reaction between PR₃ and In(1)b is more favored. Fig. 4 illustrates a schematic picture for HOMO/LUMO orbitals of CO₂, benzyne, and In(1)b in p/mCN (as a sample of electron-withdrawing groups) and p/mNHMe (as a sample of electron-releasing groups).

Fig. 4 HOMO/LUMO band gaps of CO₂, benzyne, and In(1)b in p/mCN and p/mNHMe derivatives.

In respect to the reported experimental fact,²⁴ benzyne can react with CO₂ in the absence of triphenylphosphine. On the other hand, In1(b) is more stable (−25.84 kcal mol⁻¹) than all derivatives of In1(a). It can be concluded that forming zwitterion In1(a), which accomplishes charge developing, is not a possible approach in the reaction. Thus, C–P bond formation in Path a is not an appropriate phenomenon which is an interesting question. We decided to study the phenomenon by EDA analyses. For this purpose, some triphenylphosphine derivatives, having electron-withdrawing and electron-releasing groups are chosen and relaxed scans of the C–P bond between benzyne and PPh₃ are performed, perturbing the optimized C–P bond by 0.1 Å in each step for 10 times. Then the geometry of each step (10 steps for each derivative) was extracted as an xyz format file for sobEDA analyses. Fig. 5 depicts the performed analyses for p/mCN and p/mNHMe substituents also Fig. 2(S)† depicts the analyses for the others.

Fig. 5 The EDA analyses for triphenylphosphine (H) and the derivatives p/mCN (as a sample of electron-withdrawing groups), p/mNHMe (as a sample of electron-realizing groups) in the first step of mechanism a.

3.2 EDA and deformation density analyses

Fig. 5 shows that ΔE_rep and ΔE_orb are the first and the second effective factors on ΔE_int. Also, the growth of ΔE_rep in CN, as an electron-withdrawing substitute, is more considerable than NHMe as an electron-releasing one. Generally, it can be concluded that the C–P bond formation causes more instability due to growing repulsion interactions. This negative effect is more remarkable in electron-withdrawing groups such as CN and NO₂. When it comes to ΔE_orb, as a reducing factor on ΔE_int, the electron-releasing groups induce more decline during the C–P bond development. However, ΔE_orb cannot compensate for the remarkable increase in ΔE_rep. Another interesting observation can be concluded by the compression of meta and para positions. Fig. 5 and 2(S)† depict that rising ΔE_rep values in meta positions are more significant. These reveal that the substituents in the meta position are an effective hindrance in C–P bond evolution. The comparison of diMe, Furyl, Thiopene, and Cy (Fig. 2(S)†) reveals that the effects of diMe, Furyl, and Cy are similar to electron-releasing groups, whereas Thiopen plays as a potent electron-withdrawing group. Also, as a notable observation among all of the studied PPh₃ derivatives, by C–P developing the growth of ΔE_rep in Furyl and Thiopen are the minimum and maximum values, respectively. In conclusion, regarding the discussed studies, can be deduced that the progress of the reaction via Path A is not favored, because the increase in ΔE_rep value during C–P development is a limiting factor, so other favorable factors such as ΔE_orb and ΔE_els cannot offset the unfavorable nature of ΔE_int.

The first step in mechanism b is In1(b) formation via a [2 + 2] cycloaddition reaction in which the reaction goes through a pericyclic transition state. Fig. 6 depicts the EDA study for mechanism b which reveals that both ΔE_rep and ΔE_orb have greater growth than the first step in mechanism a. Also, ΔE_els values have a more significant decline than the corresponding values in In1(a) formation in mechanism a. The sum of ΔE_els and ΔE_orb has a compensatory effect in ΔE_rep which decreases ΔE_int more remarkably than that one in mechanism a. Thus, a decrease in the ΔE_int values, as a favored factor to progress the reaction during C–C bond formation in mechanism b, is more meaningful than C–P bond development in mechanism a. It can be considered as a good justification to progress the reaction via mechanism b and not mechanism a. Fig. 7 shows the optimized TS1(b), corresponding deformation density maps and the calculated values. It is clear that Δρ^Pauli value has a more significant role in the Δρ^Total.

Fig. 6 The EDA analyses for the first step of Path b.

Fig. 7 The optimized structure of TS1(b) and the corresponding deformation density map components (isovalue = 0.005).

Step2(b) is accomplished via a nucleophilic attack of the phosphorous atom of PPh₃ to the carbon atom of In1(b) resulting in the formation of the final product. Fig. 8 shows the calculated transition states for this step. The step progression was also investigated via the EDA method. Generally, two factors can be considered effective factors in ΔE_int of the intermediate formation in which one factor decreases and another factor grows the ΔE_int. Table 3 (S)† shows the calculated EDA parameters for all the studied derivatives. The parameters such as ΔE_c, ΔE_elec, and ΔE_orb, as attraction parameters, decline ΔE_int values, however, ΔE_rep has an inverse effect which increases ΔE_int values. On the other hand, ΔE_rep values are more powerful than each attraction parameter, individually. It can be considered all of the attraction parameters as ΔE_att = ΔE_c + ΔE_elec + ΔE_orb which against ΔE_rep, has a stabilizing effect on ΔE_int values. Fig. 9 and 3(S)† depict the calculated EDA analyses for H, p/mCN and p/mNHMe. Concerning the parameters of H, it is clear that ΔE_int of p/mCN and p/mNHMe are more positive and more negative values, respectively. On the other hand, p/mNHMe possesses lower ΔE_att than the corresponding values of p/mCN. Notably, ΔE_rep for the p/mNHMe has larger values than H and p/mCN. Fig. 3(S)† depicts the same trend for other derivatives. As a result, electron-releasing groups cause an increase in ΔE_rep and a decrease in ΔE_att values, simultaneously. However, the role of electron-withdrawing groups such as CN, NO₂, CF₃, etc in altering the ΔE_rep and ΔE_att values is not as significant as that of electron-donating groups such as NHMe, OMe, Me, etc.

Fig. 8 The calculated transition states for the Step2(b).

Fig. 9 The EDA analyses H, p/mCN, and p/mNHMe in the Step2(b).

When it comes to the comparison of meta and para positions, it is obvious that the changes in ΔE_rep values for meta and para positions for electron-withdrawing and electron-releasing groups are different. Indeed, in electron-withdrawing groups changing the position from para to meta leads to a growth in ΔE_rep values. However, the change for electron-releasing groups has inverse effects on the values. In the case of ΔE_att values the trend is vice versa of ΔE_rep. As a matter of fact, changing para to meta positions causes a decline and a growth in the ΔE_att values in electron-withdrawing and electron-releasing groups, respectively.

Table 1 shows the calculated distortion (ΔE_dist) and interaction (ΔE_int) parameters for the TS2(b). Fig. 10 illustrates the plots of ΔE^≠ against the ΔE_int and ΔE_dist values in two separate categories for the meta and para positions of substitutions. These plots show a possible correlation between ΔE^≠ and ΔE_int or ΔE_dist values. Indeed, the meta-position plot includes an acceptable correlation between ΔE^≠ and ΔE_dist values (R² = 0.82). However, the para-position one depicts a remarkable correlation between ΔE^≠ and ΔE_int values (R² = 0.93). Thus, the observed ΔE^≠ values corresponding to para substituents are affected by ΔE_int, while meta position groups have their effects on ΔE^≠ values by influence on ΔE_dist. It may be related to the direct electronic effects of para substituents through the direct resonance effect which causes a higher role of ΔE_int values for this position in the C–P bond developing in TS2(b). However, in mete substituents, ΔE_dist is the determining factor because the direct resonance effect has diminished. Another fact is that ΔE_dist values in meta-position groups generally possess higher values than para position.

Table 1 ΔE^≠, ΔE_int and ΔE_dist values for the TS2(b) (all in kcal mol⁻¹)

	R	ΔE^≠	ΔEint	ΔEdist		ΔE^≠	ΔEint	ΔEdist
Para-position	NO2	27.86	−22.96	50.82	Meta-position	36.92	−26.29	63.21
	CN	26.74	−24.73	51.47		26.44	−26.6	53.04
	CF3	25.18	−27.08	52.26		39.87	−28.33	68.20
	F	20.96	−29.49	50.45		23.70	−28.54	52.24
	H	19.90	−30.3	50.20		19.90	−30.3	50.20
	Cl	22.71	−30.89	53.60		24.31	−28.64	52.95
	Me	18.71	−36.07	54.78		19.41	−34.94	54.35
	OMe	16.96	−39.03	55.99		20.55	−33	53.55
	NHMe	12.14	−44.72	56.86		17.96	−35.99	53.95
	Thiophene	23.85	−25.63	49.48		—	—	—
	Furyl	29.05	−22.36	51.41		—	—	—
	Cy	9.91	−42.45	52.36		—	—	—
	diMe	18.14	−36.7	54.84		—	—	—

---

Fig. 10 The plots of ΔE^≠ against the ΔE_int and ΔE_dist values in the meta and para positions.

Fig. 11 and 4(S)† depict the deformation density components for H, p/mCN and p/mNHMe derivatives in TS2(b). In both p/mCN and p/mNHMe the Δρ^Pauli is a higher value than Δρ^Orb which has more effects on the Δρ^Total. On the other hand, it can be found that the Δρ^Total and Δρ^Orb in p/mNHMe groups possess higher and lower values than the p/mCN ones, respectively. Thus, with respect to H, it can be concluded that electron-releasing groups vary all the components of the deformation density more remarkably than electron-withdrawing groups. When it comes to the comparison of para and meta positions, again considering the H, it is conceivable that para derivatives result in more variation for the components in both electron-releasing and electron-withdrawing groups. Indeed, in relation to the NHMe, Δρ^Orb value from −219.65 kcal mol⁻¹ in H has been reached to −229.55 kcal mol⁻¹ (ΔΔρ^Orb = −9.90 kcal mol⁻¹) in meta position, but for the para position the value reaches −242.65 kcal mol⁻¹ (ΔΔρ^Orb = −23 kcal mol⁻¹). Also, the Δρ^Pauli value in H changes from 325.12 kcal mol⁻¹ to 336.94 kcal mol⁻¹ in meta position (ΔΔρ^Pauli = 11.82 kcal mol⁻¹) and to 352.42 kcal mol⁻¹ in para position (ΔΔρ^Pauli = 27.30 kcal mol⁻¹). Similar to NHMe, for CN group variations in meta position include ΔΔρ^Orb = 8.5 kcal mol⁻¹ and, ΔΔρ^Pauli = −9.13 kcal mol⁻¹, while in the para position, the variations are ΔΔρ^Orb = 14.18 kcal mol⁻¹ and, ΔΔρ^Pauli = −15.81 kcal mol⁻¹. The same trend can be deduced in other derivatives.

Fig. 11 The deformation density components for H, p/mCN and p/mNHMe derivatives in TS2(b).

Investigation of possible correlation between deformation density components and ΔE^≠ values of TS2(b) in Fig. 12 reveals that a good correlation is observable in ΔE^≠ vs. Δρ^Pauli (R² = 0.93) and Δρ^Orb (R² = 0.94), in para positions. However, the correlation with the Δρ^Total (R² = 0.65) is not acceptable. In the meta position, all of the correlations are not satisfactory. Indeed, it seems that the substituted derivatives in para positions enforce their electronic effects as Δρ^Pauli and Δρ^Orbterms more remarkably than in the meta position.

Fig. 12 The plots of ΔE^≠ against Δρ^Pauli, Δρ^Orb and Δρ^Total in the para and meta positions.

4. Conclusion

This study provides an in-depth computational investigation of the reaction mechanisms involving triphenylphosphine (PPh₃) derivatives, benzyne, and CO₂, presenting new insights into the electronic and steric factors influencing reactivity and selectivity. Two competing mechanisms were investigated: Path a, involving direct C–P bond formation, and Path b, which progresses via a [2 + 2] cycloaddition. The results indicate that Path b is energetically favored, supported by more favorable interaction energies, and the enhanced stability of intermediates. Key findings from Energy Decomposition Analysis (EDA) and deformation density studies reveal C–P bond formation in Path a is hindered by substantial repulsive interactions (ΔE_rep), especially for electron-withdrawing substituents like CN and NO₂. The rise in ΔE_rep outweighs stabilizing factors such as orbital (ΔE_orb) and electrostatic (ΔE_els) interactions, making this pathway unfavorable. On the other hand, the formation of In1(b) through a [2 + 2] cycloaddition shows compensatory effects between attractive and repulsive forces, leading to reduced ΔE_int values. EDA analysis of TS2(b) highlights the competing influences of interaction and distortion energies. Para-substituents demonstrate a strong correlation between activation energy (ΔE^≠) and interaction energy (ΔE_int). In contrast, meta-substituents predominantly influence ΔE^≠ through distortion energy (ΔE_dist). Electron-releasing groups like NHMe enhance deformation density components such as Pauli deformation (Δρ^Pauli) and orbital deformation (Δρ^Orb) further stabilizing TS2(b). Thus, Deformation density analysis emphasizes the critical role of Δρ^Pauli and Δρ^Orb in determining transition state stability. Para-substituents exhibit stronger electronic effects, as shown by higher Δρ^Pauli and Δρ^Orb contributions compared to meta-substituents. This reinforces the importance of substituent design in modulating reactivity. Generally, this study identifies Path b as the preferred mechanism for PPh₃ derivatives, supported by both thermodynamic and kinetic favorability. The findings highlight the utility of EDA and deformation density methodologies in dissecting reaction pathways, providing mechanistic aspects that guide the rational design of PPh₃ derivatives for optimized reactivity and selectivity in CO₂ transformation reactions.

Data availability

The data underlying this study are available in the published article and its ESI.†

Conflicts of interest

The authors declare no conflict of interest.

Acknowledgements

The Research Council of the Academic Center for Education, Culture, and Research (ACECR) is gratefully acknowledged for the financial support of this project.

References

1. S. Portela and I. Fernández, J. Org. Chem., 2022, 87, 9307–9315.
2. Z. Boughlala, C. Fonseca Guerra and F. M. Bickelhaupt, J. Phys. Chem. A, 2019, 123, 9137–9148.
3. M. Baranac-Stojanović, J. Aleksić and M. Stojanović, RSC Adv., 2015, 5, 22980–22995.
4. P. Umadevi and L. Senthilkumar, RSC Adv., 2016, 6, 38919–38930.
5. S. M. N. V. T. Gorantla, M. Francis, S. Roy and K. C. Mondal, RSC Adv., 2021, 11, 6586–6603.
6. F. M. Bickelhaupt and K. N. Houk, Angew. Chem., Int. Ed., 2017, 56, 10070–10086.
7. Y. García-Rodeja, M. Sola and I. Fernandez, J. Org. Chem., 2017, 82, 754–758.
8. T. Lu and Q. Chen, J. Phys. Chem. A, 2023, 127, 7023–7035.
9. X.-F. W. Yanying Zhao and W. Li, J. Phys. Chem. A, 2024, 128, 3777–3783.
10. H. A. Rodríguez, D. A. Cruz, J. I. Padrón and I. Fernández, J. Org. Chem., 2023, 88, 11102–11110.
11. P. Mondal, N. Mandal, A. K. Pal and A. Datta, J. Org. Chem., 2024, 89, 11371–11379.
12. D. A. O'Sullivan and W. Lepkowski, Chem. Eng. News, 1990, 68, 42–61.
13. L. Merzoud, A. Saal, R. Moussaoui, O. Ouamerali, C. Morell and H. Chermette, Phys. Chem. Chem. Phys., 2018, 20, 16102–16116.
14. I. Fernández, Phys. Chem. Chem. Phys., 2014, 16, 7662–7671.
15. R. Chakraborty, J. J. Talbot, H. Shen, Y. Yabuuchi, K. M. Carsch, H. Z. H. Jiang, H. Furukawa, J. R. Long and M. Head-Gordon, Phys. Chem. Chem. Phys., 2024, 26, 6490–6511.
16. G. A. Ardizzoia and S. Brenna, Phys. Chem. Chem. Phys., 2017, 19, 5971–5978.
17. H. Sabet-Sarvestani, M. Izadyar and H. Eshghi, Energy, 2017, 134, 493–503.
18. H. Sabet-sarvestani, M. Izadyar, H. Eshghi and N. Noroozi-shad, Theoretical Approaches to CO 2 Transformations, Springer International Publishing, 2022.
19. K. J. Kron, S. J. Gomez, Y. Mao, R. J. Cave and S. Mallikarjun Sharada, J. Phys. Chem. A, 2020, 124, 5359–5368.
20. Z. Moussa, Z. M. A. Judeh and S. A. Ahmed, RSC Adv., 2019, 9, 35217–35272.
21. H. Sabet-Sarvestani, M. Izadyar, H. Eshghi and N. Norozi-Shad, Phys. Chem. Chem. Phys., 2019, 22, 223–237.
22. H. Sabet-Sarvestani, M. Izadyar, H. Eshghi and N. Noroozi-Shad, Energy, 2018, 145, 329–337.
23. H. Sabet-Sarvestani, H. Eshghi and M. Izadyar, RSC Adv., 2017, 7, 1701–1710.
24. P. Xie, S. Yang, Y. Guo, Z. Cai, B. Dai and L. He, J. Org. Chem., 2020, 85, 8872–8880.
25. P. Verma and D. G. Truhlar, Trends Chem., 2020, 2, 302–318.
26. M. Walker, A. J. A. Harvey, A. Sen and C. E. H. Dessent, J. Phys. Chem. A, 2013, 117, 12590–12600.
27. Y. Zhao and D. G. Truhlar, Theor. Chem. Acc., 2008, 120, 215–241.
28. F. Weigend and R. Ahlrichs, Phys. Chem. Chem. Phys., 2005, 7, 3297–3305.
29. G. I. M. Frisch, G. Trucks, H. Schlegel, G. Scuseria, M. Robb, J. Cheeseman, G. Scalmani, V. Barone, B. Mennucci and G. Petersson, Gaussian, 2009.
30. K. Fukui, J. Phys. Chem., 1970, 74, 4161.
31. M. Stahn, S. Ehlert and S. Grimme, J. Phys. Chem. A, 2023, 127, 7036–7043.
32. M. Rakhshanipour, H. Sabet-Sarvestani and H. Eshghi, Struct. Chem., 2020, 31, 585–598.
33. W. Z. Pengcheng Liu, J. Han, Y. Chen, H. Yu and X. Zhou, J. Phys. Chem. A, 2024, 128, 3007–3014.
34. A. Chauhan, H. S. Karnamkkott, S. M. N. V. T. Gorantla and K. C. Mondal, ACS Omega, 2022, 7, 31577–31590.
35. S. M. N. V. T. Gorantla and K. C. Mondal, ACS Omega, 2021, 6, 17798–17810.
36. S. Fakhraee and S. M. Azami, J. Chem. Phys., 2009, 130, 1–6.
37. L. Zhao, M. von Hopffgarten, D. M. Andrada and G. Frenking, Wiley Interdiscip. Rev.:Comput. Mol. Sci., 2018, 8, 1–37.
38. F. Ferraro, C. A. Barboza and E. Osorio, Chem. Phys. Lett., 2023, 811, 140240.
39. J. I. M.-A. Samir Kenouche, Nassima Bachir, ChemPhysChem, 2023, 24, e202200488.
40. O. A. Stasyuk, H. Szatylowicz, T. M. Krygowski and C. Fonseca Guerra, Phys. Chem. Chem. Phys., 2016, 18, 11624–11633.
41. T. Lu and F. Chen, J. Comput. Chem., 2012, 33, 580–592.

---

Footnote

† Electronic supplementary information (ESI) available: Supplementary tables and figures and Cartesian coordinates of the molecules have been provided in two separate files. See DOI: https://doi.org/10.1039/d5ra00804b

---