A Quantitative Molecular Orbital Perspective of the Chalcogen Bond.

We have quantum chemically analyzed the structure and stability of archetypal chalcogen-bonded model complexes D2 Ch⋅⋅⋅A- (Ch = O, S, Se, Te; D, A = F, Cl, Br) using relativistic density functional theory at ZORA-M06/QZ4P. Our purpose is twofold: (i) to compute accurate trends in chalcogen-bond strength based on a set of consistent data; and (ii) to rationalize these trends in terms of detailed analyses of the bonding mechanism based on quantitative Kohn-Sham molecular orbital (KS-MO) theory in combination with a canonical energy decomposition analysis (EDA). At odds with the commonly accepted view of chalcogen bonding as a predominantly electrostatic phenomenon, we find that chalcogen bonds, just as hydrogen and halogen bonds, have a significant covalent character stemming from strong HOMO-LUMO interactions. Besides providing significantly to the bond strength, these orbital interactions are also manifested by the structural distortions they induce as well as the associated charge transfer from A- to D2 Ch.

Introduction

The chalcogen‐bond (ChB) is the net‐attractive intermolecular interaction, often referred to as noncovalent interaction, between a Lewis‐basic chalcogen‐bond acceptor A and a Lewis‐acidic chalcogen‐bond donor D2Ch featuring a chalcogen (group 16) atom Ch to which A binds. Nearly 40 years ago, the first systematic study appeared of the chalcogen bond in which S⋅⋅⋅Y (e. g., Y = S, O, F, Cl, or Br) nonbonded atomic contacts were investigated. Early studies generally characterized chalcogen bonds as being predominantly electrostatic in nature. Later on, the significance of charge transfer from the occupied orbital of a Lewis base into an empty σ*‐type orbital of a chalcogen molecule controlling the chalcogen bond strength was recognized. Chalcogen‐bonding has since found applications in various fields of chemistry, including, supramolecular, biochemistry, spectroscopy and catalysis.

In this study, we have computationally analyzed a range of chalcogen‐bonded D2Ch⋅⋅⋅A− complexes (Ch = O, S, Se, Te; D, A = F, Cl, Br; see Scheme 1), using relativistic density functional theory (DFT) at ZORA‐M06/QZ4P. One purpose of our work is to provide a set of consistent structural and energy data from which reliable trends can be inferred for a wide range of model systems. The primary objective is to achieve a detailed understanding of the nature of chalcogen bonds by studying the associated electronic structure and bonding mechanism and compare them with the better‐known halogen bonds and hydrogen bonds.

Scheme 1

Chalcogen‐bonded D2Ch⋅⋅⋅A− model complexes (Chv=vO, S, Se, Te; D, A = F, Cl, Br).

To this end, we first explore how the geometries and energies of our model complexes D2Ch⋅⋅⋅A− vary as the chalcogen atom (Ch), or the chalcogen bond accepting Lewis base (A−) are varied. To understand the origin of the computed trends, activation strain analyses are performed on the formation of the chalcogen‐bond complexes. As part of these analyses, the interaction energy and the underlying bonding mechanism are furthermore examined in the context of quantitative Kohn‐Sham molecular orbital (MO) theory in combination with an energy decomposition analysis (EDA).[ , ] Our systematic and detailed analyses along the entire reaction profile for each of the chalcogen‐bond complexation reactions provide in‐depth insights. In particular, they demonstrate that chalcogen bonds are not at all purely electrostatic phenomena but are, to a substantial extent, covalent in nature.

Theoretical Methods

Computational Details

All calculations were carried out using the Amsterdam Density Functional (ADF) 2017.103 program. The equilibrium geometries and energies of chalcogen‐bonded complexes were computed at DFT level using the meta‐hybrid functional M06. In addition, a large uncontracted relativistically optimized QZ4P Slater type orbitals (STOs) basis set containing diffuse functions was used. The QZ4P all‐electron basis set, no frozen‐core approximation, is of quadruple‐ζ quality for all atoms and has been augmented with the following sets of polarization and diffuse functions: two 3 d and two 4 f on oxygen and fluorine, three 3 d and two 4 f on sulfur and chlorine, two 4 d and three 4 f on selenium and bromine, one 5 d and three 4 f on tellurium and iodine. The molecular density was fitted by the systematically improvable Zlm fitting scheme. The scalar relativistic effects were accounted for by using the zeroth‐order regular approximation (ZORA) Hamiltonian. It has been shown that these computational settings give accurate bond lengths and energies.

Analysis of the Bonding Mechanism

Insight into the bonding mechanism is obtained through activation strain analyses of the various chalcogen bond formation reactions. These complexation reactions are computationally modeled by increasing the distance between A− and the Ch atom of the D2Ch fragment, allowing the system to geometrically relax at each point. The D2Ch⋅⋅⋅A− distance is increased from the actual bond length value in the chalcogen‐bonded complex (r Ch⋅⋅⋅A −) to a value of 4.300 Å. Thus, each analysis starts from an optimized D2Ch⋅⋅⋅A− complex, which is then transformed to the D2Ch molecule and a halide at a relatively large distance.

These complexation reactions are analyzed using the activation strain model. The activation strain model of chemical reactivity is a fragment‐based approach to understand the energy profile of a chemical process in terms of the original reactants. Thus, the potential energy surface ΔE(ζ) is decomposed along the reaction coordinate ζ (or just at one point along ζ) into the strain energy ΔE strain(ζ), which is associated with the geometrical deformation of the individual reactants as the process takes place, plus the actual interaction energy ΔE int(ζ) between the deformed reactants [Eq. 1].

In the equilibrium geometry, that is, for ζ = ζeq, this yields an expression for the bond energy ΔE(ζeq) = ΔE strain(ζeq) + ΔE int(ζeq). The PyFrag program was used to facilitate the analyses along the reaction coordinate ζ of the bond formation processes. The interaction energy ΔE int(ζ) between the deformed reactants is further analyzed in the conceptual framework provided by the quantitative Kohn–Sham MO model. To this end, it is decomposed in three physically meaningful terms [Eq. (2)] using a quantitative energy decomposition scheme developed by Ziegler and Rauk.

The usually attractive term ΔV elstat corresponds to the classical Coulomb interaction between the unperturbed charge distributions of the deformed reactants and has four components [Eq. (3)]: i) the electrostatic repulsion between the electron densities of fragments 1 and 2, ΔV elstat,; ii) the electrostatic attraction between the nucleus of fragment 1 and the electron density of fragment 2, ΔV elstat,; iii) the electrostatic attraction between the electron density of fragment 1 and the nucleus of fragment 2, ΔV elstat,; and iv) the electrostatic repulsion between the nuclei of fragments 1 and 2, ΔV elstat,. The Pauli repulsion energy (ΔE Pauli) comprises the destabilizing interactions between occupied orbitals of the reactants and is responsible for steric repulsion. The orbital‐interaction energy (ΔE oi) accounts for charge transfer, that is, the interaction between occupied orbitals of one fragment with unoccupied orbitals of the other fragment, including the interactions of the highest occupied and lowest unoccupied MOs (HOMO−LUMO), and polarization, that is, empty–occupied orbital mixing on one fragment, due to the presence of another fragment.

The electron density distribution is analyzed using the Voronoi deformation density (VDD) method for computing atomic charges. The VDD atomic charge on atom X in a molecule (Q X VDD) is computed as the (numerical) integral of the deformation density in the volume of the Voronoi cell of atom X [Eq. (4)]. The Voronoi cell of atom X is defined as the compartment of space bounded by the bond midplanes on and perpendicular to all bond axes between nucleus X and its neighboring nuclei.

Here, the deformation density is the difference between ρ(r), i. e., the electron density of the overall molecule or complex, and ρ promolecule(r) = ΣY ρ Y(r), i. e., the superposition of spherical average‐of‐configuration atomic densities ρ Y(r) of each atom Y in the fictitious promolecule without chemical interactions, in which all atoms are considered neutral. The interpretation of the VDD charge Q Ch VDD is rather straightforward and transparent: instead of measuring the amount of charge associated with a particular atom Ch, Q Ch VDD directly monitors how much charge flows out of (Q Ch VDD > 0) or into (Q Ch VDD < 0) the Voronoi cell of atom Ch due to chemical interactions.

The VDD scheme can also be used to directly compute how much charge flows into or out of an atomic Voronoi cell X in an overall complex (e. g., [D2Ch⋅⋅⋅A]−) relative to two (poly)atomic molecular fragments (e. g., D2Ch and A−), instead of spherical atoms, as shown in [Eq. 5].

ΔQ X VDD is a measure of how the atomic charge of atom X changes due to the bonding between the fragment. In this work, [Eq. (5)] is used to compute the flow of electrons from the halide A− to the chalcogen‐bond donating molecule D2Ch (see ΔQ D2Ch VDD in Table 1).

Table 1

Activation strain analyses (in kcal mol−1) of a representative set of D2Ch⋅⋅⋅A− at the equilibrium geometries (in Å, deg.)[a]

D2Ch⋅⋅⋅A−	ΔE	ΔE strain	ΔE int	ΔQ D2Ch VDD	r Ch⋅⋅⋅A	Δr D1−Ch	Δr D2−Ch	Θ1	Θ2	ΔΘ1
F2O⋅⋅⋅F−	−21.9	28.3	−50.2	−0.37	1.784	0.408	−0.002	97.6	97.6	−6.1
F2O⋅⋅⋅Cl−	−9.9	28.1	−37.9	−0.35	2.183	0.402	0.007	98.8	101.3	−4.9
F2O⋅⋅⋅Br−	−11.5	45.0	−56.5	−0.48	2.113	0.582	0.023	98.6	103.4	−5.1
Cl2O⋅⋅⋅F−	−16.0	24.3	−40.3	−0.41	1.838	0.500	−0.032	104.4	99.5	−7.9
Cl2O⋅⋅⋅Cl−	−6.5	24.5	−31.0	−0.41	2.172	0.491	−0.012	105.5	105.5	−6.6
Cl2O⋅⋅⋅Br−	−11.0	46.5	−57.5	−0.62	1.966	0.878	0.009	105.7	110.8	−6.6
Br2O⋅⋅⋅F−	−12.9	4.5	−17.4	−0.26	2.162	0.153	0.007	106.1	86.6	−8.0
Br2O⋅⋅⋅Cl−	−6.0	1.5	−7.6	−0.19	2.673	0.084	0.021	111.9	94.0	−2.2
Br2O⋅⋅⋅Br−	−6.2	20.6	−26.7	−0.44	2.243	0.425	0.026	108.2	108.2	−5.9
F2S⋅⋅⋅F−	−50.1	16.2	−66.3	−0.35	1.813	0.227	0.045	87.0	87.0	−11.2
F2S⋅⋅⋅Cl−	−23.1	8.5	−31.7	−0.21	2.452	0.149	0.031	89.3	88.2	−8.9
F2S⋅⋅⋅Br−	−19.8	7.1	−26.9	−0.19	2.647	0.135	0.028	90.1	88.7	−8.1
Cl2S⋅⋅⋅F−	−53.1	24.5	−77.5	−0.47	1.748	0.461	0.035	92.1	93.0	−11.5
Cl2S⋅⋅⋅Cl−	−25.3	14.8	−40.1	−0.34	2.339	0.323	0.032	94.4	94.4	−9.2
Cl2S⋅⋅⋅Br−	−22.4	14.4	−36.9	−0.34	2.506	0.322	0.030	94.9	95.1	−8.7
Br2S⋅⋅⋅F−	−51.7	22.2	−73.9	−0.51	1.743	0.487	0.028	92.8	93.9	−11.1
Br2S⋅⋅⋅Cl−	−24.7	12.2	−36.8	−0.36	2.346	0.320	0.027	95.6	95.3	−8.3
Br2S⋅⋅⋅Br−	−22.3	12.3	−34.6	−0.36	2.507	0.327	0.027	96.0	96.0	−7.9
F2Te⋅⋅⋅F−	−72.4	7.9	−80.3	−0.32	2.054	0.162	0.038	84.2	84.2	−9.7
F2Te⋅⋅⋅Cl−	−42.5	5.9	−48.5	−0.24	2.608	0.134	0.035	85.2	86.3	−8.7
F2Te⋅⋅⋅Br−	−38.1	5.6	−43.7	−0.24	2.777	0.130	0.034	85.2	86.6	−8.7
Cl2Te⋅⋅⋅F−	−73.3	10.8	−84.1	−0.39	2.039	0.285	0.047	90.0	88.5	−8.2
Cl2Te⋅⋅⋅Cl−	−43.0	8.6	−51.6	−0.30	2.582	0.249	0.045	91.4	91.4	−6.8
Cl2Te⋅⋅⋅Br−	−38.6	8.4	−47.1	−0.30	2.745	0.246	0.045	91.4	91.8	−6.8
Br2Te⋅⋅⋅F−	−72.0	9.9	−81.9	−0.41	2.040	0.227	0.045	90.8	88.5	−8.3
Br2Te⋅⋅⋅Cl−	−42.0	7.6	−49.6	−0.32	2.582	0.461	0.035	92.5	92.0	−6.6
Br2Te⋅⋅⋅Br−	−37.7	7.6	−45.3	−0.32	2.751	0.487	0.028	92.8	92.8	−6.3

[a] Computed at ZORA‐M06/QZ4P. For full set of data, see Tables S1–S2 in the Supporting Information.

Results and Discussion

Chalcogen Bond Strength and Structure

The results of our ZORA‐M06/QZ4P calculations are shown in Table 1 for a representative selection of oxygen‐, sulfur‐, and tellurium‐bonded model complexes D2Ch⋅⋅⋅A−, covering D, A = F, Cl, and Br (the complete dataset for all model systems is provided in Tables S1–S2). In the first place, we note that all model reactions are associated with single‐well potential energy surfaces (PES), that is, there is no energy barrier separating the reactants from their resulting complex. In the cases where D ≠ A, C s symmetric complexes with D1−Ch bond lengths longer than D2−Ch and with bond angles Θ1 ≠ Θ2 are formed. For the cases where D = A, C 2v symmetric complexes with equal bond distances r Ch−D1 = r Ch⋅⋅⋅A and equal bond angles Θ1 = Θ2 are formed (see Table 1).

In general, chalcogen bonds become stronger on descending group 16 in the periodic table, in agreement with previous ab initio results.[ , ] The heavier D2Ch⋅⋅⋅A− chalcogen bonds (i. e., Ch = S, Se, and Te) become weaker and longer as the accepting halide (A−) varies from F− to Br−. In the case of the tellurium‐bonded complexes D2Te⋅⋅⋅A−, for example, the ΔE weakens from around −73 kcal mol−1 for A− = F− to around −38 kcal mol−1 for A− = Br− (see Table 1). However, the oxygen chalcogen bonds D2O⋅⋅⋅A− display a more complex dependency of ΔE upon variation of the accepting halide A−. From A− = F− to Cl−, the oxygen‐bond strength still weakens, similar to the situation for the heavier chalcogen bonds. However, thereafter, from A− = Cl− to Br−, the oxygen‐bond strength does not weaken but instead becomes stronger. This is most clearly seen in the series constituted by the complex F2O⋅⋅⋅A− between an oxygen molecule and a halide ion. Here, ΔE for the oxygen bond strength varies along A− = F−, Cl−, and Br− with values of −21.9, −9.9, and −11.5 kcal mol−1, respectively (see Table 1).

When the substituent D is varied from D = F to D = Br, the heavier chalcogen‐bond strength (i. e. Ch = S, Se, and Te) changes only slightly (see Table 1 and S1 in the Supporting Information). For example, along the series from F2Te⋅⋅⋅F− to Br2Te⋅⋅⋅F−, the tellurium bond strength varies only from a ΔE value of −72.4 to −72.0 kcal mol−1, the tellurium‐bond distance r Ch⋅⋅⋅A decreases in value from 2.054 to 2.040 Å, and the stretch Δr D1−Ch increases in value from 0.162 to 0.227 Å. An exception to this is again the oxygen bond D2O⋅⋅⋅A−, which becomes weaker and longer as D is varied from F to Br (see Table 1). For example, along the series from F2O⋅⋅⋅F− to Br2O⋅⋅⋅F−, the oxygen‐bond strength weakens from a ΔE value of −21.9 to only −12.9 kcal mol−1, the oxygen‐bond distance r Ch⋅⋅⋅A increases in value from 1.784 to 2.162 Å, and the stretch Δr D1−Ch decreases in value from 0.408 to 0.153 Å.

Bond Analyses with Variation of Ch

The strengthening of chalcogen bonds D2Ch⋅⋅⋅A−, as Ch varies along O, S, Se, and Te, with no change in the donating atom (D) and the accepting halide (A−), is related to the increasing electronegativity difference across the D−Ch bonds as Ch descends in the periodic table, which is translated into two main effects. Firstly, this causes the Ch atom to become increasingly positive along O, S, Se, and Te (see VDD atomic charges in Table 2), resulting in a greater electrostatic attraction. Secondly, this causes, among other effects that will be explained later, the σ* D−Ch antibonding 4a′ acceptor orbital to have higher amplitude on Ch (see Figure 1), resulting in stronger HOMO−LUMO orbital interactions.

Table 2

Bond lengths (in Å), bond angle (in deg.), VDD charge (in a.u.), orbital energies (in eV) and the homolytic bond dissociation energy without ZPE (in kcal mol−1) of isolated D2Ch fragments.[a]

D2Ch	r D−Ch	Θ1	Q Ch VDD	ϵ(1a1)	ϵ(3a1)	ϵ(1b1)	ϵ(2b1)	BDED−Ch [b]
F2O	1.376	103.7	0.09	−16.6	−2.0	−15.8	−9.6	38.0
Cl2O	1.681	112.3	−0.06	−13.8	−3.3	−13.2	−8.1	34.9
Br2O	1.818	114.1	−0.13	−12.6	−3.3	−12.0	−7.7	33.7
F2S	1.586	98.2	0.18	−15.0	−0.7	−14.5	−7.1	92.2
Cl2S	2.016	103.6	0.12	−12.1	−2.5	−11.4	−7.1	64.0
Br2S	2.180	103.9	0.04	−11.3	−2.7	−10.5	−7.0	51.4
F2Se	1.730	96.3	0.28	−13.9	−1.6	−13.4	−7.2	87.6
Cl2Se	2.155	101.1	0.21	−11.6	−2.7	−10.8	−7.0	63.2
Br2Se	2.312	101.7	0.13	−10.8	−2.8	−10.1	−6.9	50.6
F2Te	1.892	93.9	0.31	−13.0	−2.0	−12.7	−6.6	92.8
Cl2Te	2.333	98.2	0.29	−11.0	−2.6	−10.3	−6.6	69.3
Br2Te	2.492	99.1	0.22	−10.3	−2.7	−9.6	−6.4	55.7

[a] Computed at ZORA‐M06/QZ4P; [b] Energy for the reaction D2Ch→DCh.+D..

Figure 1

Schematic molecular orbital diagram for a) isolated D2Ch fragments at C 2v symmetry (blue: a1; green: a2; red: b1; black: b2) and b) D2Ch⋅⋅⋅A− complexes. The first column in (b) refers to the isolated D2Ch fragment and the second column refers to the D2Ch fragment deformed to its C s symmetric geometry in the complex (blue: a′; red: a′′), in which one D−Ch bond has been elongated. See Figure S1 in the Supporting Information for the 3D isosurfaces of the orbitals.

The trend in bond strength ΔE is mainly determined by the interaction energy ΔE int. For example, from F2O⋅⋅⋅F− to F2Te⋅⋅⋅F−, ΔE is strengthened from −21.9 to −72.4 kcal mol−1 while ΔE int is strengthened from −50.2 to −80.3 kcal mol−1 (see Table 1). The trend in ΔE is further enhanced by the strain energy (ΔE strain), which becomes less destabilizing from Ch = O to Te. However, the differences are smaller than the differences in ΔE int. For example, from F2O⋅⋅⋅F− to F2Te⋅⋅⋅F−, ΔE strain is weakened by 20.4 kcal mol−1 (from 28.3 to 7.9 kcal mol−1; see Table 1), while ΔE int becomes 30.1 kcal mol−1 more stable.

To understand the origin of these trends, we have carried out activation strain analyses along the entire reaction coordinate ζ, projected onto the stretch in D1−Ch bond, Δr Ch−D1, that occurs as the chalcogen‐bond accepting A− atom approaches the D2Ch molecule (see Theoretical Methods section). The resulting activation strain diagrams (ASD) including EDA terms of the interaction are shown for a representative example series, namely, F2O⋅⋅⋅F− to F2Te⋅⋅⋅F−, in Figure 2 (for the complete dataset, see Table S2). Again, the trend in bond energy ΔE is mainly determined by ΔE int(ζ), which strengthens when going from Ch = O to Te (Figure 2, left). On the other hand, the ΔE strain(ζ) curves almost coincide. However, the strain curves reach a final point at ζeq, that is, the equilibrium geometry of the complex; and here the strain energy ΔE strain(ζeq) becomes more destabilizing from Ch = Te to O. Note that the trend in strain energies at the equilibrium geometries along the series of F2Ch⋅⋅⋅F− complexes (see Table 1) arises mainly from changes in the steepness of the interaction curves, not from the relatively minor variation in the strain curves (see Figure 2). Thus, as the F2Ch⋅⋅⋅F− interaction gets weaker along Ch = Te, Se, S and O, the interaction curve becomes shallower and the balance between strain and interaction curve, i. e., the stationary point of the complex, occurs at longer and longer F−Ch distances and, consequently, more destabilizing ΔE strain(ζeq) (see Table 1).

Figure 2

Activation strain (left panel) and energy decomposition (right panel) analyses of F2Ch⋅⋅⋅F− chalcogen‐bonded complexes (green, Ch = O; black, Ch = S; blue, Ch = Se; red, Ch = Te).

To understand the trends in ΔE int(ζ), we further decomposed the ΔE int into the individual energy components (Figure 2, right). The electrostatic energy ΔV elstat(ζ) is the least stabilizing for Ch = O and then strengthens along S, Se, and Te. This can be understood by the increasing differences in electronegativity across the D−Ch bonds when going from O to Te, resulting in a larger positive charge on Ch. For example, the VDD atomic charge on Ch in F2O, F2S, F2Se, and F2Te amounts to +0.09, +0.18, +0.28, and +0.31 a.u., respectively, and becomes even more positive as the D1−Ch bond elongates (see Figure 3a). Nevertheless, our analyses reveal that the chalcogen bonding mechanism is absolutely not purely electrostatic but instead has a relatively large covalent component (ΔE oi), stemming mainly from the HOMO−LUMO interaction between the occupied halide npy atomic orbital (AO) and the σ* D−Ch antibonding 4a′ acceptor orbital (see Figure 1). The associated charge transfer from A− to D2Ch is reflected by the ΔQ D2Ch VDD, which is negative, i. e., D2Ch gains charge from A− upon complexation, for all D2Ch⋅⋅⋅A− complexes (see Table 1). For example, ΔQ D2Ch VDD is −0.37 a.u. for F2O⋅⋅⋅F− and −0.32 a.u. for F2Te⋅⋅⋅F−. The HOMO−LUMO charge transfer nature of the chalcogen bond is also clearly reflected by the associated deformation density. This is illustrated by the 3D plots of the deformation densities associated with chalcogen‐bond formation in F2S⋅⋅⋅F− and F2Te⋅⋅⋅F− (see Figure 4). As can be seen, there is charge depletion on the Lewis base F− (and in between the Ch⋅⋅⋅F− bond due to the Pauli repulsion ) and charge accumulation on D2Ch. Note the 3D shape of the regions of charge depletion and accumulation: they reflect the shape of the 2p‐type lone pair from which the F− Lewis base donates and the σ* D−Ch antibonding 4a′ acceptor orbital on D2Ch into which this charge is donated, respectively, in the HOMO−LUMO interaction. For the chalcogen bonded complexes, the orbital interaction term ranges from 37 % for F2Te⋅⋅⋅F− to as much as 76 % for Br2O⋅⋅⋅F− of the total bonding interactions (ΔE oi + ΔV elstat; see Table S2 in the Supporting Information). As can be seen in our energy decomposition diagram, the orbital interaction curves ΔE oi(ζ) become more stabilizing from Ch = O to Te (Figure 2, right). The stronger orbital interaction for the heavier chalcogens is the result of the larger LUMO−HOMO overlap (i. e. ⟨4a′|npy⟩; see Figure 1 for the MO diagram which shows the npy orbital of A− pointing towards the D1−Ch bond of the D2Ch fragment) as Ch becomes more electropositive. For example, in the Cl2Ch⋅⋅⋅Cl− series, ⟨4a′|npy⟩ increases from 0.12 to 0.20 to 0.22 to 0.24 along Ch = O, S, Se, and Te in the equilibrium geometry (see Table S2 in the Supporting Information). The larger percent contribution of the covalent component on oxygen bonds is simply because the electrostatic attraction is relatively weak, caused by the smaller positive charge on O (see Table 2).

Figure 3

Energy of the 4a′ orbital (in eV) and the VDD charge on Ch atom (in a.u.) in the neutral fragment D2Ch projected onto a) the F1−Ch bond stretch (green, Ch = O; black, Ch = S; blue, Ch = Se; red, Ch = Te) and b) the D1−Te bond stretch (black, D = F; blue, D = Cl; red, D = Br).

Figure 4

Deformation density (Δρ(r)=ρ [D2Ch⋅⋅⋅A‐](r)−ρ D2Ch(r)−ρ A‐(r); red = depletion; blue = accumulation) plot (a and c) and HOMO−LUMO interaction (b and d) for a representative series of D2Ch⋅⋅⋅A− chalcogen bonds.

Whereas ΔE oi(ζ) becomes more stabilizing from Ch = O to Te, it becomes comparable in magnitude for all chalcogens in the equilibrium geometry ΔE oi(ζeq). This is a consequence of the fact that the F1−Ch bond expansion becomes more pronounced when going from Ch = Te to O. The increasing F1−Ch bond expansion causes the σ* D−Ch antibonding 4a′ acceptor orbital (see Figure 1a) to drop further in energy for lighter chalcogens, resulting in a smaller HOMO−LUMO gap and hence more stabilizing donor‐acceptor interactions. This effect can be observed in Figure 3a, which shows the energies of the σ* F−Ch antibonding 4a′ acceptor orbitals along the reaction coordinate. For Ch = S, Se, and Te, the energy of the σ* F−Ch antibonding 4a′ acceptor orbital converges to an energy value of −3.8 eV as the chalcogen bond is formed. For Ch = O, on the other hand, the σ* F−O antibonding 4a′ acceptor orbital energy quickly drops to a value of −6.4 eV, because the overlap between the F and O AOs is more sensitive to the D−Ch distance than for the more diffuse AOs of heavier Ch (see Figure S2 in the supporting information). However, it is counteracted by the orbital overlap between the σ* D−O antibonding 4a′ acceptor orbital and the npy donor orbital, which is significantly worse for Ch = O than for other chalcogen systems (see Table S2 in the Supporting Information).

Bond Analyses with Variation of A−

Our analyses show that the weakening of heavier chalcogen bonds D2Ch⋅⋅⋅A− (Ch = S, Se, Te), as the accepting group varies from A− = F− to Br−, is directly related to the concomitant reduction in electron‐donating capacity of the np‐type HOMO and thus the Lewis basicity of the A− halide. We recall the chalcogen bonds display both an electrostatic component (ΔV elstat) and a covalent component (ΔE oi). The latter stems mainly from the HOMO−LUMO interaction between the occupied halide np atomic orbital (AO) and the σ* D−Ch antibonding 4a′ acceptor orbital (see Figure 1). Both ΔV elstat and ΔE oi are weakened as the halide HOMO becomes more diffuse and effectively lower in energy from A− = F− to Br− (see Table S2). Consequently, the interaction energy (ΔE int) and, thus, the net chalcogen‐bond strength ΔE becomes less stabilizing along A− = F− to Br− (see Table 1 and Table S1 in the Supporting Information). This is very similar to what was found for hydrogen bonds DH⋅⋅⋅A− and heavier halogen bonds DX⋅⋅⋅A− (X = Cl, Br, I).

The key to understanding why oxygen bonds D2O⋅⋅⋅A− show a more complex, partially opposite trend (i. e., the expected weakening from A− = F− to Cl− but thereafter a strengthening along A− = Cl− to Br−) is contained in the counteracting effects evolving from D−O bond stretching induced in the triatomic D2O molecule as it interacts with the halide A−. Interestingly, activation strain analyses reveal, again, that interaction energies recover the original trend in total energies, that is, ΔE int(ζ) weakens from A− = F− to Br−. This can be seen in Figure 5 which shows the activation strain and energy decomposition diagrams along the reaction coordinate ζ projected onto the stretch Δr D1−Ch for two representative series. Each diagram in Figure 5 refers to one particular F2O or F2Te molecule forming a chalcogen bonding with A− = F−, Cl−, and Br−. The ΔE strain curves within each subgraph coincide because they refer to the same D−Ch bond in the same triatomic molecule being stretched as the complexation reaction progresses. Consequently, the trend A− = F− to Br− in the total F2O⋅⋅⋅A− and F2Te⋅⋅⋅A− energy profiles ΔE in each subgraph is directly determined by the trend in the corresponding ΔE int curves.

Figure 5

Activation strain (left panel) and energy decomposition (right panel) analyses of a) F2O⋅⋅⋅A− and b) F2Te⋅⋅⋅A− (black, A− = F−; blue, A− = Cl−; red, A− = Br−).

The reason why the oxygen bonds D2O⋅⋅⋅A− do not experience a weakening in ΔE int from A− = F− to Br−, as all other chalcogen bonds, is promoted by a combination of factors: i) a weak D−O bond that is easily stretched; ii) a strong interaction with an approaching halide A−; and iii) a σ* D−Ch antibonding 4a′ acceptor orbital that drops in energy, more quickly than for other D1−Ch bonds due to a more sensitive overlap between the D1 and O AOs, as the D1−O bond elongates (see Figure 3 and Figure S2 in the supporting information). The latter generates a stronger driving force for D1−Ch stretching in D2Ch⋅⋅⋅A− because this deformation enhances the orbital interactions and thus ΔE int. Note that, for D2O⋅⋅⋅A−, ΔE oi is the strongest bonding component and that the ΔE oi(ζ) curves directly reflect the electron‐donating capacity of the np‐type HOMO of the A− halides, that is, the ΔE oi curves become more stabilizing from A− = Br− to F− (see Figure 5). Indeed, D1−Ch stretching is most pronounced if this bond in the neutral fragment is weaker, that is, for the weaker chalcogen bonds (e. g., ca. 38 kcal mol−1 for F−O, ca. 35 kcal mol−1 for Cl−O and ca. 34 kcal mol−1 for Br−O; see Table 2). In this case, it is able to affect the trend in overall bond strength ΔE. The D1−O stretching in oxygen‐bonded complexes is most pronounced in the Cl2O⋅⋅⋅A− series, along which the Cl1−O stretch Δr D1−Ch varies between 0.5 and 0.9 Å, but it is already relevant in the F2O⋅⋅⋅A− series in which the F1−O stretch Δr Ch−D1 varies between 0.4 and 0.6 Å from A− = F− to Br− (see Table 1).

We conclude that, in general, chalcogen bonds D2Ch⋅⋅⋅A− become weaker along A− = F− to Br− because the larger radii and lower np AO energies of the halides lead to weaker electrostatic attraction and weaker orbital interactions. The trend in D2O⋅⋅⋅A− oxygen bond strength is partially inverted, that is, ΔE becomes more stabilizing along A− = Cl− and Br− because of a subtler interplay of factors. Notably, a significant stretching of the relatively weak D−O bonds in the D2O⋅⋅⋅A− equilibrium structures lowers the σ* D−O antibonding 4a′ acceptor orbital and thus amplifies the donor‐acceptor orbital interactions.

Bond Analyses with Variation of D

The strength of the heavier chalcogen bonds D2Ch⋅⋅⋅A− varies little when going from D = F to Br because the Cl−Ch and Br−Ch bonds are significantly weaker than the F−Ch. This allows the Cl−Ch and Br−Ch bonds to stretch to a higher extent and, therefore, to have more stabilizing electrostatic attraction and orbital interactions. For the oxygen bonds D2O⋅⋅⋅A−, the bond energy is weakened along the same variation because the D−Ch bond strength are all comparable (see Table 2). In both cases, the trend in bond strength ΔE is determined by the interaction energy ΔE int. For example, from F2O⋅⋅⋅F− to Br2O⋅⋅⋅F−, ΔE int is weakened from −50.2 to −17.4 kcal mol−1, respectively, whereas from F2Te⋅⋅⋅F− to Br2Te⋅⋅⋅F−, the bond energy only changes from −80.3 to −81.9 kcal mol−1, respectively (see Table 1). The strain energy (ΔE strain) is not negligible, but it does not offset the trend set by ΔE int. Our activation strain analyses explain the above differences between oxygen and heavier chalcogen bonds (see Figure 6).

Figure 6

Activation strain (left panel) and energy decomposition (right panel) analyses of a) D2O⋅⋅⋅F− and b) D2Te⋅⋅⋅F− (black, D = F; blue, D = Cl; red, D = Br).

Starting with some general observations, we find that for oxygen, as well as heavier chalcogen bonds, the ΔE strain curves are most unfavorable when D = F and gradually become less destabilizing as the donating atom is varied along D = F, Cl, and Br (see Figure 6). Furthermore, for all D2Ch⋅⋅⋅A− complexes, the ΔE int curves become less stabilizing along D = F, Cl, and Br. The resulting energy profiles of D2Ch⋅⋅⋅A− depend on the balance between both ΔE strain and ΔE int, but the interaction energy curves already show a very similar trend to ΔE.

The slope and shape of the ΔE strain curves is of course directly related to the D1−Ch bond strength of the neutral fragment, which in general becomes stronger as the polarity across the D−Ch bond increases (see Table 2). From F2Ch to Br2Ch, where Ch is S, Se or Te, the halogen‐chalcogen bond strength decreases significantly from a value of ca. 93 to 50 kcal mol−1 (Table 2). The corresponding halogen−oxygen bonds are all much weaker, and variations in the homolytic bond dissociation energy (BDE) are also much smaller. From F2O to Br2O, the bond strength decreases from 38.0 to 33.7 kcal mol−1. Thus, for the heavier chalcogen‐bonded complexes, where Ch is S, Se, or Te, the ΔE strain curves show a pronounced reduction in slope from F2Ch to Br2Ch, which, in the corresponding chalcogen‐bonded complexes F2Ch⋅⋅⋅A− to Br2Ch⋅⋅⋅A−, translates into an increasing stretch Δr D1−Ch of the neutral fragment. As the stretch Δr D1−Ch becomes larger from equilibrium structures F2Ch⋅⋅⋅A− to Br2Ch⋅⋅⋅A−, the ΔE int curves have been able to descend further, to lower, more stabilizing energies. This stabilization is, of course, related to the ΔV elstat and ΔE oi. Note that the electrostatic attraction and orbital interaction curves become less stabilizing along D = F, Cl, and Br, but turn out to have comparable strength in the equilibrium structures, because the D1−Ch bonds have been increasingly stretched in the latter, that is, in Cl2Ch⋅⋅⋅A− and Br2Ch⋅⋅⋅A−. The bonding components ΔV elstat and ΔE oi are the most stabilizing for D = F because of the larger difference in electronegativity across the D−Ch bonds (vide supra). However, the ΔE oi is able to further stabilize for D = Cl and Br because, in the equilibrium structure of the chalcogen‐bonded complexes, the Cl−Ch and Br−Ch bonds expand to a higher extent, resulting in a stronger stabilization of their σ* D−Ch antibonding 4a′ acceptor orbitals (see Figure 3b). Furthermore, the VDD atomic charge on Ch becomes increasingly more positive as the D1−Ch bond expands, which translates into more stabilizing ΔV elstat for D = Cl and Br in the equilibrium geometry. The final result is, thus, a comparable stability among heavier chalcogen bonds D2Ch⋅⋅⋅A− complexes when the substituent D is varied from F to Br.

Chalcogen Bonds Versus Halogen and Hydrogen Bonds

Our analyses highlight that chalcogen bonds, halogen bonds, and hydrogen bonds are all similar in nature. Each of these bonds in our set of model systems has a significant covalent component in addition to electrostatic attraction, and can range in strength roughly between −6 and −70 kcal mol−1 (see Figure 7). Chalcogen bonds and halogen bonds have a larger range in polarities in D−Ch and D−X than in D−H bonds and are in general stronger than hydrogen bonds because of more stabilizing orbital interactions (see Table S3 for bond energies ΔE of a representative series of XB and HB ΔE computed at ZORA‐M06/QZ4P). However, chalcogen bonds and halogen bonds also have more destabilizing Pauli repulsion because the lone‐pair HOMO of the Lewis base overlaps with more closed shells, in particular, with the σ D−Ch bonding 3a′ and 2a′ FMOs or σ D−X bonding FMO with a higher amplitude on Ch and X, respectively, than the amplitude of σ D−H bonding FMO has on H (see Figure 7; see also Ref. 9). Our analyses provide a unified picture for chalcogen bonds, halogen bonds, and hydrogen bonds based on quantitative Kohn‐Sham molecular orbital theory, which proves that these intermolecular interactions cannot be described by a pure and simple electrostatic model.

Figure 7

Generic molecular orbital diagrams for a) D2Ch⋅⋅⋅A− chalcogen bonds, b) DX⋅⋅⋅A− halogen bonds, and c) DH⋅⋅⋅A− hydrogen bonds.

Conclusions

Chalcogen bonds in D2Ch⋅⋅⋅A− range between 6 and 73 kcal mol−1 in strength, becoming stronger as the chalcogen atom becomes more electropositive, along Ch = O, S, Se and Te, and also as the halide becomes a stronger Lewis base, along A− = Br−, Cl− and F−. The trend upon variation of the substituent along D = F, Cl, Br is less pronounced, as are all trends for the relatively weak oxygen bonds. This follows from our bonding analyses based on relativistic density functional theory (DFT) calculations at ZORA‐M06/QZ4P.

Our activation strain and quantitative Kohn‐Sham MO bonding analyses reveal that the chalcogen bonds in D2Ch⋅⋅⋅A− are similar in nature to halogen bonds in DX⋅⋅⋅A− and hydrogen bonds in DH⋅⋅⋅A− (Ch = O, S, Se, Te; D, X, A = F, Cl, Br). Chalcogen bonds are far from being solely electrostatic phenomena. Similar to halogen and hydrogen bonds, chalcogen bonds have a sizeable covalent component, ranging up to 80 % of the bonding components (ΔV elstat + ΔE oi), stemming from HOMO−LUMO interactions between the np‐type lone pair on the bond accepting fragment A− and the LUMO with strong D−Ch σ* anti‐bonding character on the bond donating fragment D2Ch.

Chalcogen bonds become stronger for heavier Ch because of the greater difference in electronegativity across the D−Ch bonds, causing: i) the σ* D−Ch antibonding 4a′ acceptor orbital to have higher amplitude on Ch, enhancing HOMO−LUMO orbital interactions; and ii) the Ch to become more positively charged, resulting in greater electrostatic attraction when descending in group 16 of the periodic table. The chalcogen bonds also become stronger for lighter A− because the electron‐donating capacity of the np‐type HOMO (i. e. Lewis basicity) of the halides increases ascending group 17 in the periodic table. The trends for oxygen bonds, as well as along various D, are less pronounced because of counteracting effects or small values in bond strength.

Conflict of interest

The authors declare no conflict of interest.

As a service to our authors and readers, this journal provides supporting information supplied by the authors. Such materials are peer reviewed and may be re‐organized for online delivery, but are not copy‐edited or typeset. Technical support issues arising from supporting information (other than missing files) should be addressed to the authors.

Supplementary

Click here for additional data file.