Topochemical reactions are a kind of solid-state reaction that depend strongly on the stacking of the reactant molecules.1,2 Initiated by heating, irradiation or compression, such reactions can be used to produce a variety of interesting molecules and polymeric materials with desirable optoelectronic and mechanical properties.3–5 They are also highly atom-economic, with excellent regio- and stereo-selectivity, and they sometimes result in single-crystal-to-single-crystal transformation.6,7 Since the first systematic exploration of the [2 + 2] reaction of substituted cinnamic acids by Schmidt in 1964,2,8 more and more topochemical reactions have been developed, including [4 + 4], [4 + 2], and [3 + 2] cycloaddition reactions, as well as 1,4- or 1,6-addition reactions of diynes, dienes, trienes, and triynes.1,7,9–12 All these investigations suggest that suitable stacking of the molecules and short spacing between the reactive groups are key to topochemical reactions, and some topochemical postulates have been proposed to provide guidance for the design of the reactions.2 For instance, the distance between C=C double bonds in [2 + 2] cycloaddition reactions should be in the range of 3.5–4.2 Å, while the ideal distance for the topochemical polymerization of diacetylene is d_(C1⋯C4ʹ) < 4 Å, with stacking spacing ≈5 Å.7,8 However, most of the proposed distances have been determined from the original and ground-state crystal structures, rather than the structure at the moment of the photo- or thermal reactions, and therefore they are not closely related to the intrinsic reaction criteria, nor are they really predictive.13,14

Recently, pressure-induced polymerization (PIP) has attracted widespread attention, given that the application of pressure can continuously adjust the intermolecular distance and drive molecules to spontaneously react, where a critical structure before the reaction is clearly defined.15 A uniform empirical shortest intermolecular C⋯C distance at threshold, d_(C⋯C)min, can be derived from PIP reactions at room temperature (RT), suggesting that it is directly related to the nature of the functional groups involved.16 As summarized in Table I, although the reaction pressure varies, alkynes always react when d_(C⋯C)min reaches ∼3 Å. For example, acetylene dicarboxylic acid polymerizes to form crystalline percarboxylpolyacetylene at d_(C⋯C)min ≈ 3.1 Å, which is almost the same value as for acetylene.5,17 By contrast, aromatics react at ∼2.8 Å owing to their relative inertness,20,21 and alkynylarenes have a longer distance owing to their concerted bonding.22,23

This empirical threshold d_(C⋯C)min is an excellent parameter to predict topochemical reactions under external pressure. However, its physical and chemical basis is still unknown. A quantitative understanding of this threshold is critical for designing rational syntheses via topochemical reactions under pressure and is even important for solid-state heat/irradiation-initiated reactions. An investigation of the PIP of benzene crystals established a link between the maximum fluctuation of d_C⋯C between rigid molecules on the one hand and the temperature T and the phonon frequency ν (expressed as wavenumber) of the lowest optical translational mode at the Γ-point of the Brillouin zone on the other, with the square of the molecular displacement x being given byx2=kBTM(2πcν)2,(1)where k_B is Boltzmann’s constant, M is the molecular mass, and c is the speed of light.20 By subtracting the thermal fluctuations from the experimentally determined distance between rigid molecules under different pressure–temperature conditions, it was found that the distance converged to ∼2.6 Å as a trigger threshold. This work revealed the important role of phonons in PIP, but the underlying physical mechanism remained obscure, and the complex reaction pathways of benzene cannot be fully understood from this model.

Acetylene is the simplest alkyne, and its polymerization is a textbook addition reaction that initiated the study of conductive polymers, making it a perfect model to study the origin of the threshold distance and the reaction selectivity.24 It was reported that acetylene polymerized to trans-polyacetylene at 4 GPa and RT,25 or to cis-polyacetylene at 12 GPa and 77 K.26 In the most recent crystallographic investigation, acetylene was found to polymerize at 5.7 GPa and RT via a topochemical route, forming cis-polyacetylene and subsequently nano-graphane.17In situ neutron diffraction and molecular dynamics (MD) simulations disclosed that the bonding started at d_(C⋯C)min = 3.1 Å. In the present paper, we adopt a computational approach to probe the threshold distance of acetylene, decompose the reaction process of PIP, reveal the underlying physics, and uncover the origin of the RT experimental threshold distance of 3.1 Å. These quantitative results should also help in understanding other topochemical addition bonding process.

To investigate the structural variations under applied pressure, we first optimized the crystal structure of acetylene at 0 K and external pressure. Starting with the experimental critical structure at RT and 5.7 GPa reported in Ref. 17, stress-equilibrated structures were obtained under different pressures. In our calculation, the optimized structure at 0 GPa had the closest lattice parameters and d_(C⋯C)min to the experimental critical structure (experimental 5.7 GPa; for a comparison, see Table S1 in the supplementary material), and so we used the optimized result at 0 GPa to represent the experimental crystal structure at 5.7 GPa. We investigated the pressure range from 0 to 30 GPa and found that acetylene polymerized above 28.8 GPa. Representative structures and geometric parameters are displayed in Fig. 1. The acetylene molecules lie on the b–c plane. With increasing pressure, the cell parameters shrink with a decreasing rate until the acetylene is about to polymerize. a is the most compressed, followed by c, and b is the least. For the two closest intermolecular distances, d1 in the a–c plane (potentially leading to the cis-polymer) decreases faster than d2 before reaction. In addition, the angle ∠C–C–c (∠1) continues to increase, while ∠H–H–c (∠2) does not do so before the reaction occurs, which results in bending of the H–C≡C–H molecule [smaller bond angle ∠3, marked in Fig. 1(a)]. The equilibrated C≡C bond lengthens by about 1% just before reaction, indicating a weakening.

At 28.9 GPa, acetylene polymerized during optimization. We then concluded that the d_(C⋯C)min (d1 = 2.3 Å) at 0 K and 28.8 GPa was the intrinsic threshold distance of this PIP. We performed an MD simulation at 28.8 GPa and a limited temperature (50 K) to visualize the bonding process, in which the chain polymer was first generated along the a + c/a − c direction (Fig. S1, supplementary material), consistent with Ref. 17. For the crystal structures obtained during optimization and MD simulation, we studied the relationship between intermolecular electronic interaction and d_C⋯C by single-point calculations. The Mayer bond order of the intermolecular C⋯C bond and the value of the electron localization function (ELF) at the midpoint of the bond, which respectively represent the sharing of electron pairs and the electronic localization between these two carbon atoms, are plotted vs bond length in Fig. 2. Both curves have turning points near 2.3 Å, suggesting enhanced intermolecular electron interaction at that distance. This weakens the repulsion between adjacent molecules and corresponds to the tipping point of bonding. Thus, we can conclude that intermolecular bonding in the acetylene crystal commences when the 2.3 Å intrinsic threshold distance is reached.

To probe the dynamic properties of the acetylene crystal and find the possible reactions involved, we first calculated the phonon spectrum of the structure at 0 K and threshold pressure for reference. On forcing the last stress-equilibrated optimized structure before reaction (28.8 GPa) to a little higher pressure (29 GPa) with limited optimization steps, an imaginary frequency appeared at the Γ-point and T-point (0.5, 0, −0.5), as shown in Fig. 3(a). The imaginary mode at the Γ-point (−76 cm⁻¹) belongs to an optical branch originating from the lattice vibrational mode (Fig. S2, supplementary material), which is −44 cm⁻¹ at the T-point. Another three vibrational modes with imaginary frequencies are also found at −44 cm⁻¹ at the T-point, corresponding to the acoustic branches. This tells us that both the vibrations with the T-point symmetry and those with the Γ-point symmetry (symmetry of vibrations, or collective behaviors of molecular motion) should be highly focused.

In RT experiments, the static external pressure changes the equilibrated configuration, while it is always the thermal displacement of atoms leading to a decrease in d_C⋯C that initiates the bonding. Therefore, we calculated the phonon spectrum of the structure at 0 GPa [experimental 5.7 GPa, Fig. 3(c)] and analyzed the vibrational modes in detail to help demonstrate the effect of temperature. Figure S3 (supplementary material) presents details of 48 normal vibrational modes at the Γ-point, which are closely related to molecular motions. Modes with wavenumber below 300 cm⁻¹ are lattice vibrational modes (modes 1–20 in Fig. S3). Modes with wavenumbers in the range 600–900 cm⁻¹ (modes 21–36) are bending modes. Stretching modes for C–C bonds (modes 37–40) and C–H bonds (modes 41–48) are located around 2000 and 3100 cm⁻¹, respectively. Among these vibrations, lattice vibrational modes, corresponding to translations or rotations of the four molecules, are the major contributors to the fluctuations in d_C⋯C. According to the analysis in Sec. II, we focus on the decrease of d1. Let us take the translations along a as an example (modes 1–4 in Fig. S3). In modes 1 (acoustic branch) and 2, the adjacent carbon atoms from neighboring molecules move in the same direction with d_C⋯C unchanged, whereas in modes 3 and 4, they move in opposite directions to each other and compress d_C⋯C. As a result, about half of the lattice vibrational modes, specifically modes 3, 4, 11, 12, 15, 16, 19, and 20 (group 1) at the Γ-point, play the major roles in decreasing the value of d_C⋯C. By contrast, modes 1, 2, 9, 10, 13, 14, 17, and 18 (group 2) have little effect. For the vibrational modes with T-point symmetry, it is still the lattice vibrational modes that cause d_C⋯C to decrease, which will be illustrated later. Thus, the vibrations closely related to intermolecular bonding have been determined.

In solids, the collective behaviors of atoms correspond to phonon modes. For a system of N atoms in the solid state, there are 3N vibrational degrees of freedom, corresponding to 3N phonon branches. For simplicity, a classical physical model is used here when dealing with the total thermal energy of the system, i.e., this is given by 3Nk_BT at a relatively high temperature T according to the Dulong–Petit law.

When we focus on a certain Q-point in the Brillouin zone (or a certain wave vector) from the viewpoint of classical physics, the 3N corresponding vibrational modes at this point, with frequencies given by the phonon dispersion relations, are mutually independent. In this situation, the average thermal energy of each mode is the same, giving k_BT in total. Thus, under the harmonic approximation, the maximum displacement of N atoms in the jth phonon mode at this Q-point can be represented as∑iNatoms12mi‖xi,j,max⇀‖2=kBT(2πcνj)2,(2)

where xi,j,max⇀ is the maximum displacement of the ith atom determined by the jth phonon mode, ν_j represents the frequency of this mode (expressed as wavenumber), and m_i is the mass of the ith atom. Equation (2) is an enhancement of Eq. (1) and provides a general illustration of how phonon modes contribute to atomic motion in solids. It also demonstrates that low-frequency vibrations, lattice modes for example, contribute more to the decrease in d_C⋯C than high-frequency vibrations such as stretching modes. With Eq. (2), we can estimate the instantaneous shortest d_C⋯C in the acetylene system at RT.

In the limiting case, the maximum displacement of a carbon atom isxC,max⇀=∑jxC,j,max⇀.(3)

The instantaneous shortest d_C⋯C is obtained when all the modes contributing to the reduction in d1 are applied together. On the basis of the above analysis, both lattice vibrational modes satisfying Γ-point symmetry and those satisfying T-point symmetry can induce a significant decrease in d_C⋯C. At the Γ-point, it is the vibrational modes in group 1 that are active in this respect. However, for the vibrational modes at the T-point, there are extra phase differences in the molecular vibrations compared with the modes at the Γ-point (Fig. S4, supplementary material). Both vibrational modes in group 1 and those in group 2 at the T-point can induce a decrease in d_C⋯C, but these reductions occur in different diagonal directions (a + c/a − c) respectively. Owing to their lower frequencies, the modes in group 2 lead to larger atomic displacements of carbon atoms at the T-point; see Fig. S5 (supplementary material), which shows in greater detail part of the lattice vibrational modes from the Γ-point to the T-point in Fig. 3(c).

The maximum atomic displacements of carbon atoms induced by the modes in group 2 at the T-point were calculated according to Eq. (2) and are listed in Table II. According to Eq. (3), each of the two carbon atoms can move ∼0.4 Å to its limiting position via thermal vibrations, with the translational vibrations contributing about two-thirds of this movement and the librational vibrations about one-third. Interestingly, the movement is precisely along the C⋯C direction, as shown in Fig. 4, and so the maximum fluctuation of d_C⋯C is the sum of two opposite displacements, which is ∼0.8 Å. Thus, d_C⋯C decreased from 3.1 to 2.3 Å, the same as the aforementioned intrinsic threshold distance at 0 K. Parallel investigation of the Γ-point situation (Table S2, supplementary material, which should be compared with Table II) shows that the maximum displacement is ∼0.3 Å for one carbon atom, which is smaller than that at the T-point and is therefore less important to the PIP.

The positions of the hydrogen atoms in Fig. 4 were determined by an additional optimization with the carbon atoms frozen since these were also affected by other vibrational modes such as bending modes. We thereby obtained a transient “threshold model” structure by adding eight vibrational displacements to the experimentally determined crystal structure at RT and 5.7 GPa. In the model, the minimum instantaneous d_C⋯C reproduced the local structure of the intrinsic threshold model well at 0 K and 28.8 GPa. Thus, we could conclude that the experimentally observed threshold distance 3.1 Å is actually the sum of an intrinsic threshold of 2.3 Å and a vibrational contribution of 0.8 Å.

Finally, we performed variable cell nudged elastic band (vc-NEB) calculations to investigate the relationship between the transient “threshold model” and the transition state of PIP from the perspective of enthalpy. First, the minimum energy path (MEP) from acetylene crystal to cis-polyacetylene was searched in a 2a× 1b× 2c supercell. The energy curve is displayed in Fig. 5. In the state with the highest enthalpy (like a transition state, “TS”), the instantaneous shortest intermolecular C⋯C distance is about 1.6 Å within a dimer. The activation energy is around 2.7 eV per supercell, which is ∼134 kJ/(mol C–C bond). As polymerization proceeds, cis-polyacetylene chains are formed, accelerating the reaction of other molecules. Finally, cis-polyacetylene is formed. Based on the MEP, we constructed the pathway via the transient “threshold model,” also shown in Fig. 5 (red line). On this pathway, other structures between the initial state and the “TS” were regenerated by linear interpolation of the lattice parameters and the atomic positions. In the MEP, the unit cell changes first, and only once the intermolecular spacing has become short enough does the reaction begin (Fig. S6, supplementary material). On the new pathway, however, when 2.3 Å is reached via vibrations, intermolecular bonding is facilitated, leading to an almost flat energy curve in the bonding stage. This process is similar to the spontaneous reaction at 0 K and 28.8 GPa and demonstrates that vibrations can provide enough energy for the molecules to overcome the “TS” if they move to the vibrational maximum at RT.

An electronic interaction analysis of the transient “threshold model” also supports the similarity of the two reaction processes. The intermolecular bond order is 0.15, and the ELF value at the midpoint is 0.25, which is in accord with tipping point at 28.8 GPa in Fig. 2, implying a similar electronic interaction. Thus, we can conclude that owing to thermal vibrations, the d_(C⋯C)min under experimental conditions (3.1 Å) is dynamically compressed to the intrinsic threshold distance for PIP (2.3 Å), inducing the bonding of acetylene molecules as a trigger of the reaction. Additionally, all the structures in the two reaction processes were spin-unpolarized (the transient “threshold model” for example; see Fig. S7, supplementary material). In other words, no single electron appeared in the reaction processes. The same conclusion was reported in the benzene system in Ref. 20.

MD simulations also were carried out to provide an intuitive proof. The result of Pair distribution function (PDF) evolution are displayed in Fig. S8 (supplementary material), where a threshold for acetylene at about 2.3 Å can be observed, in excellent agreement with the conclusion above.

With the above understanding, we can now predict the influence of temperature quantitatively. Without taking changes in frequencies into consideration, the vibrational displacement is approximately proportional to the square root of temperature. Therefore, the contributions of thermal vibration are 0.9, 1.0, and 0.7 Å for 400, 500, and 200 K, respectively. The experimental threshold distances then become 3.2 Å at 400 K, 3.3 Å at 500 K, and 3.0 Å at 200 K, which are instructive values with regard to the selective synthesis of cis- or trans-polyacetylene.

In summary, our study has revealed an intrinsic threshold distance 2.3 Å for acetylene crystals at 0 K, where intermolecular C⋯C interaction is significantly enhanced. We have found that lattice vibrations are crucial for pressure-induced polymerization, calculated their contributions to the decrease in d_C⋯C at room temperature, and have applied this finding to the experimentally determined threshold structure. We have obtained a transient “threshold model” with an instantaneous d_C⋯C of 2.3 Å, which perfectly reproduces the intrinsic threshold and thus demonstrates how the reaction is initiated by vibrations. This model, which can be described as “intrinsic threshold + vibration maximum,” perfectly explains the experimental threshold of 3.1 Å for acetylene, and can also be applied on other alkyne systems. Our work has provided a quantitative understanding of the distance threshold and the effect of temperature in topochemical reactions, and the results obtained here will be useful for developing novel topochemical reactions as well as tailored syntheses of regio- and stereospecific polymeric materials.

Acknowledgment. The authors acknowledge the support of the National Natural Science Foundation of China (NSFC) (Grant Nos. 22022101, 21875006, 11704024, and 12174200). The authors also acknowledge support from the National Key Research and Development Program of China (Grant No. 2019YFA0708502) and from the Nature Science Foundation of Tianjin (Grant No. 20JCYBJC01530).