Theoretical Design of Stable Pentacoordinate Boron Compounds.

Through theoretical computations, we found that boron can form thermodynamically stable pentacoordinate compounds. Pentacoordinate boron (penta-B) is just hypercoordinate but not hypervalent because it forms only four covalent bonds, of which at least one is a multicenter bond. Being electron deficient, to be pentacoordinate, at least two of its bonding atoms should have low electronegativity. Penta-B can be formed in H k B(CH3) m (XH3) n (X = Si, Ge, Sn, and n ≥ 2) and BR5 (R = BH2NH3, AsH2, and BeH). Based on a systematic investigation of these model compounds, we designed three thermodynamically stable penta-B compounds that can potentially be synthesized by hydrogenating their tricoordinate counterparts under mild reaction conditions.

Introduction

Hypercoordination is the property of the main-group elements in a molecule having a larger than normal coordination number, typically greater than four.[1] Hypercoordination is common for the elements in period 3 and beyond.[2] However, it is difficult to form hypercoordinate compounds for the elements in period 2.[3] Over the past few decades, several research studies on hypercoordination in period 2 elements with more than three valence electrons such as carbon and nitrogen have been reported.[2d,4] However, for electron-deficient elements such as boron, there is still debate over whether they can form hypercoordinate compounds and what their bonding nature is if they exist.

The attempts on synthesizing hypercoordinate single-boron-center compounds were unsuccessful from our point of view. Since 1984, a series of the so-called pentacoordinate boron (penta-B) compounds has been synthesized, by forcing a tricoordinate boron center to form two additional bonds with Lewis-base ligands.[5,6] However, the two additional B–X (X = O, N, or Cl) bonds in these compounds (B–O: ∼2.4 Å; B–N: ∼2.5 Å; and B–Cl: ∼2.7 Å) are much longer than normal. In addition, Wiberg bond indexes (WBI)[7] of the B–X bonds are all below 0.15. Thus, they can hardly be regarded as real penta-B compounds because there are no covalent bonds formed between B and the other two additional ligands. So far, only one theoretical study mentioned five hypothetical silylboranes whose boron centers look like real pentacoordinate.[8] However, no electronic structure analyses and thermodynamic properties were provided, and it is unknown what their bonding nature is and whether they are thermodynamically stable. Because normal tricoordinate silylboranes can be synthesized and have many interesting properties,[9] we wonder if it is possible to synthesize penta-B compounds from normal tricoordinate silylboranes. Hence, in this work, we first studied the electronic and geometric properties of a hierarchy of model penta-B compounds to reveal their bonding nature. Then, we try to design several thermodynamically stable silylboranes with a penta-B center that may be synthesized by experiments under mild reaction conditions.

Results and Discussion

We first studied the electronic structure and stability of penta-B silylboranes in detail. The geometries of 10 structurally stable silylboranes, HB(CH3)(SiH3) (k = 1∼5, m = 0∼2, n = 1∼5, and k + m + n = 5), optimized using the M06-2X/aug-cc-pVTZ method,[10] are listed in Figure . To be pentacoordinate, the five bonds around B should have normal covalent bond lengths of the ordinary B–H, B–Si, and B–C bonds, which are about 1.20, 2.03, and 1.56 Å, respectively. The distances of the longest B–Si/B–H bond in each HB(CH3)(SiH3) molecule (results for the B–C bonds are not listed because they all have normal covalent bond lengths) are tabulated in Table . The B–Si/B–H bonds in 1b∼1e, 2b∼2d, and 3a are all shorter than 2.09/1.24 Å, while in 1a and 2a, they are longer than 2.20/1.30 Å. Therefore, 1b∼1e, 2b∼2d, and 3a can be regarded as penta-B compounds from a geometrical point of view.

Figure 1

Geometries of HB(CH3)(SiH3) (k = 1∼5, m = 0∼2, n = 1∼5, and k + m + n = 5).

Table 1

Longest B–Si/B–H Bond Lengths in HB(CH3)(SiH3) (k = 1∼5, m = 0∼2, n = 1∼5, and k + m + n = 5)

	bond length (Å)		bond length (Å)
1a	2.20/1.30	2a	2.27/1.35
1b	2.05/1.22	2b	2.09/1.24
1c	2.02/1.21	2c	2.05/1.23
1d	2.03/1.21	2d	2.08/–
1e	2.04/–	3a	2.06/–

WBI[7] analysis results show that the weakest B–Si/B–H bonds in 1a and 2a are just 0.38/0.62 and 0.33/0.54, respectively. They are certainly neither completely broken nor normal covalent bonds. Both geometrical and WBI data suggest that 1a and 2a are η2-complexes formed through the interaction between a σ Si–H bond orbital of SiH4 and the 2p empty orbital of the B center. On the other hand, most of the B–Si/B–H bonds in the other eight molecules can be viewed as weak covalent bonds because their WBIs are in the range of [0.47, 0.84]/[0.77, 0.93] (Table S4). Therefore, they can be regarded as penta-B compounds.

The low WBIs of the B–Si/B–H bonds of the eight molecules suggest that they are not normal single covalent bonds. Meanwhile, the sum of the WBIs of five B–X (X = H, Si, or C) bonds for the eight molecules is no more than 8, implying that no more than four covalent bonds are formed. These results suggest that B must have formed multicenter bonds with its five bonding atoms. Adaptive natural density partitioning (AdNDP) analyses[11] show that there are four, three, and two multicenter bonds around the B center in 1c, 2c, and 3a (Figure ; Figure S1 for the other five molecules), respectively. These results suggest that boron is not hypervalent,[1] and the Lewis octet rule is not violated. Boron accommodates five bonding atoms by forming multicenter bonds with some of them.

Figure 2

AdNDP multicenter orbitals of 1c, 2c, and 3a, where the pink, yellow, cyan, and white balls represent B, Si, C, and H atoms, respectively, and ON is the corresponding occupation number of the AdNDP orbital. H atoms in the silyl and methyl groups are omitted for clarity.

For formation of multicenter bonds, at least two of the bonds between boron and its ligands should not be too strong. It is critical for boron to be hypercoordinate. In addition, because boron is electron deficient, forming covalent bonds requires its bonding atoms to give electrons to boron. Therefore, at least two bonding atoms must have similar electronegativity to boron. Other than silicon, beryllium, germanium, tin, and arsenic also meet such criteria. Replacing the Si atoms by Ge/Sn in HB(CH3)(SiH3) (n ≥ 2), 16 penta-B compounds (Figure ) could be optimized. Similarly, 3 penta-B compounds, BR5 (R = BH2NH3, AsH2, and BeH), could be optimized.

Figure 3

Geometries of HB(CH3)(XH3) (X = Ge, Sn, k = 1∼5, m = 0∼2, n = 2∼5, and k + m + n = 5) and BR5 (R = BH2NH3, AsH2, and BeH).

Hypercoordination usually means instability. Table tabulates the Gibbs free-energy changes (ΔrG) and barriers (ΔG‡) of five decomposition reactions of 2b (Table S5 for the other seven) at 298.15 K. The results indeed show that these hypothetical compounds are unstable. Among them, releasing SiH4 is the easiest one with a ΔG‡ of just 2.9 kcal/mol. Releasing disilane is the second-easiest one, with a 13.8 kcal/mol ΔG‡ and a negative ΔrG. However, the positive ΔrG of the pathway to release H2 and its low ΔG‡ imply that it is possible to synthesize pentacoordinate silylboranes by hydrogenating their tricoordinate counterparts with at least two silyl groups. To obtain pentacoordinate silylboranes stable at room temperature (RT), we need to increase ΔG‡ of the lowest-energy decomposition pathway so that they are kinetically stable at RT or increase ΔrG to be positive if ΔG‡ must be low so that they are thermochemically stable.

Table 2

Gibbs Free-Energy Changes (ΔrG) and Barriers (ΔG‡) of the Five Decomposition Pathways of 2b at 298.15 K (in kcal/mol) in the Gas Phasea

decomposition products of 2b	ΔG‡	ΔrG
(SiH3)2BCH3 + H2	12.0	0.8
H2B(SiH3)CH3 + SiH4	2.9	2.3
H2BCH3 + Si2H6	13.8	–8.6
HB(SiH3)2 + CH4	20.3	–0.4
H2BSiH3 + H3SiCH3	23.1	–2.5

The results were computed by the G4//M06-2X/aug-cc-pVTZ method.

The first attempt we tried is to replace the SiH3 groups in HB(CH3)(SiH3) by more realistic SiR3 (R = methyl (Me) or phenyl (Ph)) groups. A total of five such compounds (“A” series) were designed (A1 to A5, Figure S2). Among them, A5 (pentacoordinate H2B(SiPh3)3) is a potential candidate of thermodynamically stable penta-B compounds. To design it, we make use of the π–π stack interaction to stabilize A5 and to increase ΔrG and ΔG‡ of the decomposition pathways. In addition to the π–π stack effect, another effect could be the electron-withdrawing effect of the phenyl group, which makes Si more electron deficient for forming such hypercoordinate bonding. Table tabulates the ΔGrs and ΔG‡s of three possible decomposition pathways of A5 at 298.15 K in heptane, the solvent to synthesize B(SiPh3)3.[9a] It should be noted that solvation effect is small for these reactions. The results are indeed promising because all three decomposition pathways have positive ΔrG. Considering that the pathway to HB(SiPh3)2 has a low ΔG‡ and boranes with at least one B–H bond may dimerize, we computed the ΔrG of the reaction 2A5 → (HB(SiPh3)2)2 + 2HSiPh3 and 3H2 + 2B(SiPh3)3 → 3(SiPh3)2 + B2H6. For the first reaction, we found that it has a positive ΔrG of +1.4 kcal/mol. For the other reaction, although it has a very negative ΔrG of −37.3 kcal/mol, the second step of the reaction, i.e., to release Si2Ph6 from A5, has too high a ΔG‡ (39.7 kcal/mol) to occur at RT. These results show that A5 is indeed thermodynamically stable. In addition, the B(SiPh3)3 + H2 → A5 reaction has a negative ΔrG (ΔrGH) and a moderate ΔG‡ (ΔGH‡) at RT. ΔGH‡ is so low that it is possible to synthesize A5 by hydrogenating B(SiPh3)3 even below RT.

Table 3

Gibbs Free-Energy Changes (ΔrG) and Barriers (ΔG‡) of Three Possible Decomposition Pathways of A5 at 298.15 K (in kcal/mol) in Heptane Solutiona

decomposition products of A5	ΔG‡	ΔrG
B(SiPh3)3 + H2	31.9	18.3
HB(SiPh3)2 + HSiPh3	15.6	12.3
H2BSiPh3 + Si2Ph6	39.7	9.4

The energies were computed by the M06-2X functional.

However, a potential disadvantage of synthesizing A5 is that H2B(SiPh3)3 has two low-energy conformers in fast equilibrium (Figure ). Except BLYP and B3LYP, other seven functionals (PBE, ωB97XD, M06-L, MN15-L, M06, M06-2X, and MN15) all give similar results that A5 is 1.6∼3.6 kcal/mol lower in free energy than η2-H2B(SiPh3)3 at 298.15 K (Table S9). Although an ∼2 kcal/mol free-energy lowering cannot guarantee that the experiment can surely obtain A5 other than η2-H2B(SiPh3)3 due to uncertainties in theoretical computations, our results indicate that the chance to observe A5 by hydrogenating B(SiPh3)3 is high and it is worth a try because B(SiPh3)3 had been synthesized in 1984 with a relatively easy method.[9a]

Figure 4

Gibbs free-energy profiles of the B(SiPh3)3 + H2 reaction at 298.15 K in heptane solution computed by the M06-2X functional.

Other than the “A” series, we have designed other 17 penta-B compounds (Figure S2): the “B” series (B1 to B7) containing two silyl groups, and the “C” series (C1 to C7) containing three silyl groups. Backbones were used to constrain the silyl groups and hinder the release of HSiR3 and (SiR3)2. Their stability and tricoordinate counterparts (removing two bonding H atoms on B) have been studied by searching all possible decomposition and deformation pathways (Figure S3) based on knowledge of chemical bonding and reactions. The ΔrG and ΔG‡ values of the lowest-energy pathway are named as ΔrGMin or ΔGMin‡, respectively. Promising penta-B candidates should have (1) negative ΔrGH and low ΔGH‡, i.e., they are easy to be synthesized from hydrogenating their tricoordinate counterparts (Table S10), and (2) positive ΔrGMin or high ΔGMin‡, i.e. both the pentacoordinate silylborane and its tricoordinate counterpart are stable at RT (Table S11). For ΔG‡, we set the criteria for ΔGH‡ to be better below 25 kcal/mol and ΔGMin‡ to be better above 25 kcal/mol, based on an estimation of the half-life of reaction from classical transition state theory. The calculations show that the half-life of a unimolecular reaction is about 66 h and the half-life of a bimolecular reaction is about 94 h at 298.15 K with a ΔG‡ of 25 kcal/mol. The geometries of two promising compounds, B3_Me from the “B” series and C5 from the “C” series that meet such criteria, are presented in Figure . Between them, we would recommend the synthesis of B3_Me first because it has fewer backbones and consequently can be synthesized more easily.

Figure 5

Geometries of two stable pentacoordinate silylboranes. Hydrogen atoms on the carbons are omitted for clarity.

Once these pentacoordinate silylboranes are synthesized, they can be verified by NMR spectroscopy. Table tabulates the NPA charges and the 11B NMR chemical shifts (δ(B)) of simple and recommended silylboranes. Penta-B draws electrons from silyl groups and are negatively charged. Consequently, it is much more shielded than tricoordinate boron and has very negative δ(B) values. Indeed, for penta-B silylborane compounds, δ(B) has a linear relationship with the NPA charge of the boron atom (Figure S4). On the other hand, neutral tricoordinate silylboranes have very positive δ(B) values. In addition, η2-complex and pentacoordinate conformers can be well differentiated by NMR spectroscopy because δ(B) of the pentacoordinate conformer is more negative than that of the η2-complex conformer: For H2B(SiPh3)3, δ(B) of A5 is −51.92 ppm, while that of η2-H2B(SiPh3)3 is −37.31 ppm. The same phenomenon can be observed for the hypothetical silylboranes, 1a and 2a. They are two η2-complexes and their δ(B) are just −42.04 and −23.65 ppm, respectively.

Table 4

11B NMR Chemical Shifts (δ(B)) and NPA Charges of the B Center in Selected Silylboranesd

compound	δ(B) (ppm)	NPA charge of B (a.u.)
1a	–42.04	–0.54
1b	–57.01	–1.13
1c	–67.37	–1.51
1d	–72.29	–1.72
1e	–78.11	–1.92
2a	–23.65	–0.20
2b	–48.39	–0.75
2c	–60.39	–1.14
2d	–65.33	–1.36
3a	–55.50	–0.84
B(SiPh3)3	155.15	–0.18
A5a	–51.92	–1.40
η2-H2B(SiPh3)3b	–37.31	–1.18
B3_Me (-H2)c	112.94	0.12
B3_Me	–41.63	–0.87
C5 (-H2)c	152.66	–0.39
C5	–56.71	–1.87

Pentacoordinate conformer of H2B(SiPh3)3.

η2-complex conformer of H2B(SiPh3)3.

Tricoordinate counterpart removing two bonding hydrogen atoms on B.

δ(B) was computed by a scaling method at the mPW1PW91/6-311+G(2d,p) level of theory.[12]

Conclusions

In summary, boron can be pentacoordinate by forming multicenter covalent bonds with elements, e.g., Be, B, Si, Ge, Sn, and As, having similar electronegativities to boron. We showed that penta-B is not hypervalent and does not violate the Lewis octet rule. Although hypercoordination usually implies instability, we designed three thermodynamically stable pentacoordinate silylboranes, A5, B3_Me, and C5, that may potentially be synthesized by hydrogenating their tricoordinate counterparts under mild reaction conditions. Potential usage of pentacoordinate silylboranes recommended in this study is hydrogenation catalysts or reductants.

Computational Methods

Validation of Computational Methods

First, we took pentacoordinate silylborane H3B(SiH3)2 (Figure ) as an example to test the effect of basis set and method on the geometry. Eight methods including MP2, M06-2X,[10a] MN15,[13a] ωB97XD,[13b] PBE0,[13c] TPSSH,[13d] DSD-PBEP86,[13e,13f] and PBE0DH[13g] were tested. Two basis sets were tested: (a) a large basis set in which the aug-cc-pVTZ (AVTZ) basis set[10b−10d] was used for all the atoms; (b) a smaller basis set (SBS) in which the 6-31+G(d,p) basis set[14a,14b] was used for B and its five bonding atoms, whereas the 6-31G(d,p) basis set[14c−14g] was used for the other atoms. Selected bond lengths of H3B(SiH3)2 are listed in Table . From the results in Table , it can be concluded that the method and basis set both have small effect on the geometry of H3B(SiH3)2: The standard deviation of all 16 combinations of methods and basis sets is below 0.01 Å for each B–X (X = H or Si) bond; the mean absolute deviation between the two basis sets for the eight methods is also below 0.01 Å for each B–X bond.

Figure 6

Geometry of H3B(SiH3)2 and the indices of atoms. H atoms in the silyl groups are omitted for clarity.

Table 5

Selected Bond Lengths of H3B(SiH3)2 (Unit: Å)

	B1–Si2	B1–Si3	B1–H4/H5	B1–H6
M06-2X/SBSa	2.049	2.007	1.221	1.197
M06-2X/AVTZb	2.037	2.000	1.218	1.197
MN15/SBS	2.038	1.996	1.221	1.198
MN15/AVTZ	2.025	1.984	1.217	1.194
ωB97XD/SBS	2.041	2.009	1.224	1.202
ωB97XD/AVTZ	2.031	2.001	1.222	1.203
PBE0/SBS	2.029	2.008	1.223	1.202
PBE0/AVTZ	2.024	2.002	1.222	1.203
TPSSH/SBS	2.050	2.016	1.225	1.200
TPSSH/AVTZ	2.040	2.009	1.225	1.203
MP2/SBS	2.049	2.009	1.214	1.192
MP2/AVTZ	2.034	2.004	1.218	1.196
DSDPBEP86/SBS	2.050	2.010	1.222	1.198
DSDPBEP86/AVTZ	2.041	2.005	1.225	1.200
PBE0DH/SBS	2.031	2.002	1.221	1.198
PBE0DH/AVTZ	2.022	1.996	1.222	1.200
standard deviation	0.009	0.007	0.003	0.003
MADc	0.008	0.006	0.002	0.003

SBS: 6-31+G(d,p) for B and its five bonding atoms and 6-31G(d,p) for other atoms.

AVTZ: aug-cc-pVTZ.

Mean absolute deviation between the two basis sets for the eight methods.

Second, we tested the performance of density functional theory (DFT) methods on computing relative energies. We used the M06-2X/AVTZ method to optimize the geometries of eight pent-B silylboranes, HB(CH3)(SiH3) (k = 1∼5, m = 0∼2, n = 2∼5, and k + m + n = 5). A total of 23 decomposition reactions (Table S2) for them were studied. The geometry optimizations and harmonic vibrational frequency analyses of the reactants, transition states, and products were all performed with the M06-2X/AVTZ method. The G4 method[15] was used to perform single-point energy calculations on these optimized geometries. The errors of the M06-2X method using two basis sets, AVTZ and SBS, on the energetics of the 23 reactions are summarized in Table . The results in Table indicate that the smaller basis set (SBS) systematically underestimates the reaction energies and energy barriers by about 1 kcal/mol. On the other hand, using a larger basis set, AVTZ, there is almost no systematic error because the mean error is close to 0. In addition, using AVTZ, both reaction energies and energy barriers are improved and the overall root mean square error (RMSE) over 46 relative energies is just 1.2 kcal/mol. Therefore, these results indicate that geometry optimization can be performed using an SBS, whereas a larger basis set is better to be used to further refine the energetic results. In the present study, for small systems, we used the M06-2X/AVTZ method to optimize geometry and used the G4 method to perform single-point energy calculations. For large systems where G4 calculations are prohibitively expensive, we used M06-2X/SBS to optimize geometry and a larger basis set to perform single-point energy calculations. The basis set used for single-point energy calculations is also a combined basis set: the AVTZ basis set was used for B and its five bonding atoms, and the aug-cc-pVDZ (AVDZ) basis set was used for the other atoms. This combination basis set was abbreviated as LBS. LBS is a compromise between accuracy and efficiency because for those key atoms involved in bond making and broken processes, the large AVTZ basis set was used, whereas for other “observing” atoms, a smaller AVDZ basis set was used.

Table 6

Mean Error (ME), Mean Unsigned Error (MUE), and Root Mean Square Error (RMSE) of the M06-2X Method Using the AVTZ and SBS Basis Sets on Computing All 46 Relative Energies of the Reaction Using G4 as the Standard (Unit: kcal/mol)

	totala
	M06-2X/SBSb	M06-2X/AVTZc
ME	–1.1	0.0
MUE	1.2	1.0
RMSE	1.6	1.2

The G4 results were computed at the G4//M06-2X/AVTZ level.

SBS: 6-31+G(d,p) for B and its five bonding atoms and 6-31G(d,p) for other atoms. Geometries were fully optimized using this basis set.

AVTZ: aug-cc-pVTZ. Geometries were fully optimized using this basis set.

Although G4 is very accurate, it is prohibitively expensive for large systems. We should find a cheaper method as accurate as possible. In the present study, we tested a total of 11 DFT functional methods: two GGA functionals, BLYP and PBE; two hybrid GGA functionals, B3LYP and ωB97XD; two meta-GGA functionals, M06-L and MN15-L; three hybrid meta-GGA functionals, M06, M06-2X, and MN15; and two double hybrid GGA functionals, PBE0DH and DSD-PBEP86. All energetic results were obtained by performing single-point energy calculations on geometries optimized with the M06-2X/AVTZ method. The errors of these functionals using G4 as the standard are summarized in Table . The results in Table indicate that M06-2X is the best method, which has the smallest error on both reaction energy changes and energy barriers. Therefore, for the calculations of large systems, we will use the M06-2X/SBS method to perform both geometry optimizations and vibrational frequency analyses. Then, we will use the M06-2X/LBS method to perform single-point energy calculations.

Table 7

Mean Error (ME), Mean Unsigned Error (MUE), and Root Mean Square Error (RMSE) of All 46 Relative Energies of 11 DFT Methods Using G4 as the Standard (Unit: kcal/mol)a

	BLYP	PBE	B3LYP	ωB97XD	M06-L	MN15-L	M06	M06-2X	MN15	PBE0DH	DSD-PBEP86
ME	–8.2	1.1	–6.9	–1.9	0.8	0.5	–0.8	0.0	2.1	–1.5	–8.0
MUE	8.3	2.5	7.0	2.1	2.4	1.3	1.9	1.0	2.3	2.2	8.3
RMSE	10.5	3.0	9.1	2.8	3.2	1.6	2.2	1.2	3.0	3.3	9.7

All 46 relative energies: 23 potential energy changes and 23 potential energy barriers of the reaction. Single-point energy calculations were performed on geometries optimized with the M06-2X/AVTZ method.

Computational Details

All quantum calculations were performed with the Gaussian 16 program package.[16] Wiberg bond indices,[7] NPA atomic charges, and AdNDP orbitals[11] were performed with the NBO[17] and Multiwfn programs.[18] In all the DFT calculation, a pruned (99,590) grid (using keyword “int. = ultrafine” in Gaussian 16) was used. Solvation effect was considered using the polarizable continuum solvation model[19a] with radii and nonelectrostatic terms for Truhlar and co-workers’ SMD solvation model.[19b] Except H2, for which it is a gas under standard state, Gibbs free energy of a compound in the solution (GT) is computed by the following equation:where Ee is the electronic energy computed with solvation effect considered, ZPE and ΔG0 → T are the zero-point vibrational correction and thermal correction to Gibbs free energy in the gas phase computed by the M06-2X/SBS method, and 1.9 kcal/mol is the correction from the gas-phase standard state of 1 bar to the solution standard state of 1 mol/L. The geometry used for the single-point calculation in the solution was optimized in the gas phase.

11B NMR chemical shifts were computed by a well validated method.[12] In this method, NMR calculations were computed with the mPW1PW91/6-311+G(2d,p) method in THF solution under the SMD solvation model. The 11B NMR chemical shift (δ(B), in ppm) was computed using the following scaling equation:where σ is the computed isotropic shielding constant, and the intercept and slope are 106.67 ppm and −1.1050, respectively.