Precise equilibrium structure determination of thiophene (c-C₄H₄S) by rotational spectroscopy—Structure of a five-membered heterocycle containing a third-row atom

Recently, a very satisfying level of agreement was achieved for pyrimidine between the semi-experimental equilibrium (r_e^SE) parameters and the highest practical level of theoretical calculation available—differences no more than 0.0004 Å in bond distances and no more than 0.01° in bond angles.¹ That result raises the question of whether or not such a level of agreement is possible in parallel treatments of other similarly sized molecules. While we have begun the study of similar species consisting of carbon, nitrogen, and hydrogen, we have also chosen to explore the impact of including third-row atoms in this methodology. Heavier atoms not only inherently present additional challenges regarding the computational level of theory and basis set but also may impact the ability of the theory to adequately treat the corrections to the rotational constants. The present study of thiophene was undertaken with the primary goal of attempting to answer that question for a prototypical sulfur-containing aromatic ring. Each component of this investigation represents a challenging undertaking: (i) chemical synthesis and purification of a substantial family of isotopologues of thiophene; (ii) measurement, assignment, and fitting of the high-resolution rotational spectrum of each isotopologue; (iii) computationally demanding calculations involving electronic structure, vibration–rotation coupling, and corrections for electron mass distribution; and (iv) theoretical methods for providing a “best theoretical estimate” (BTE) of the molecular structure. Integration of all of these data affords a highly precise semi-experimental structure (r_e^SE) for thiophene. The method steps beyond the traditional Kraitchman analysis of structure determination by single-atom isotopic substitution and validates the ability of modern experimental and computational methods to establish the structure of an organic heterocyclic molecule containing a third-row element.

The five-membered, aromatic ring of thiophene (c-C₄H₄S, Fig. 1, C_2v) is a fundamental molecular structural moiety in organic chemistry and is prevalent in pharmaceuticals² and organic electronic materials.^3,4 Thiophenes are naturally occurring on the Earth and have been detected on Mars.^5,6 Thiophene itself has been proposed to be present in the interstellar medium,⁷ and unsaturated carbon-chain molecules containing sulfur have been detected.^8,9 The rotational and vibrational spectroscopies of thiophene have been of interest for over a century. The first infrared absorption spectrum was reported in 1905 by Coblentz,¹⁰ and other Raman and infrared studies were reported in the 1930s and 1940s.^11–14 In 1965, a low-resolution infrared study of liquid and vapor phase thiophenes and all its deuterio isotopologues (mono- through tetra-deuterio, except [2,4-²H]-thiophene) was published, including the assignment of their fundamental vibrations.¹⁵ The vibrational assignments were later reaffirmed,¹⁶ and the rotationally resolved IR spectrum of the ν₁₃ band was reported by Pankoke et al.¹⁷ The rotationally resolved IR spectrum has subsequently been analyzed for ν₁₄, ν₈, and ν₇ (400–750 cm⁻¹)^18,19 and ν₁₃, ν₇, ν₆, ν₅, and ν₁₉ (600–1200 cm⁻¹).²⁰ These studies provided the rotational and centrifugal distortion constants for each of these vibrationally excited states, along with their precise vibrational frequencies. More recently, Raman and resonance Raman spectra of thiophenes, neat and in solution (cyclohexane, dimethyl sulfoxide, and methanol), have also been reported.^21,22 Analysis of the fundamental, combination, and overtone vibrational states using rotational spectroscopy is on-going and will be the subject of a future work.

The dipole moment of thiophene (0.55 ± 0.01 D)²³ permits its observation by rotational spectroscopy. Bak et al. provided the rotational constants for the normal isotopologue of thiophene, four deuterated isotopologues,²⁴ the ³⁴S, and the two ¹³C isotopologues based on microwave spectral data from 12 to 30 GHz.²⁵ Using Kraitchman’s equations,²⁶ they determined the first substitution structure of thiophene (r_s) from the normal and five mono-substituted isotopologues.²⁵ These early microwave data were combined with electron diffraction and liquid-crystal nuclear magnetic resonance (NMR) data to determine the r_α structure.²⁷ Later, the microwave, electron diffraction, and vibrational spectroscopy were combined with computational methods to find a mixed-method, semi-experimental equilibrium structure (r_e^SE).²⁸ The work of Bak et al.^24,25 was used to determine a subsequent semi-experimental equilibrium structure (r_e^SE) with corrections for vibration–rotation interaction constants and electron mass distributions from density functional theory (DFT) methods.^29,30 Although not utilized in the previous thiophene structure determinations, additional rotational transitions of the normal isotopologue have been reported in later works,^31–34 as well as transitions for the ³³S and ³⁴S isotopologues.³⁴ We have extended the number of isotopologues for which rotational constants are available to 26 isotopologues. Of these, 24 were included in the determination of the semi-experimental equilibrium structure (r_e^SE) presented in this work. We have systematically investigated the extent to which additional isotopologues are valuable in reducing the statistical errors of the r_e^SE structural parameters.

Rotational spectra of thiophene were collected at room temperature using three different spectrometers. The spectrum of thiophene was collected at the University of Toyama using two different instruments. The spectrum was collected between 40 and 121.5 GHz at sample pressures of 30–50 mTorr using instrumentation described previously.^35,36 Additional transitions were measured between 8 and 18 GHz with a CP-FTMW spectrometer, based on the design of McJunkins and Brown,³⁷ with a chirp span of 50–150 or 50–200 MHz, at sample pressures of 9–13 mTorr, with 10 000–50 000 molecular signal acquisitions. The broadband spectra of thiophene and its deuterium-containing isotopologues were collected at the University of Wisconsin–Madison from 130 to 230 and 235–360 GHz with a pressure of 5 mTorr using instrumentation described previously.^38,39 The spectrum of the sample containing primarily tetradeuterio thiophene was collected only from 130 to 230 GHz. Spectral segments were combined into a single broadband spectrum using Kisiel’s Assignment and Analysis of Broadband Spectra (AABS) software package.^40,41 ASFIT/ASROT were used for least-squares fits and spectral predictions,^40,41 along with PLAN and AC programs for analysis.⁴² Measured transitions for each isotopologue were least-squares fit in the I^r representation and the A and S reductions. The ground state of the normal isotopologue has also been fit in both reductions in the III^r representation. With the exception of a few new CP-FTMW measurements with 30 kHz uncertainties, all new measurements assume an experimental uncertainty of 50 kHz; all previously reported measurements use the uncertainty that was provided in either the original publication or a related publication using the same instrument.

Density functional theory calculations were performed using Gaussian 16⁴³ with the WebMO interface.⁴⁴ The optimized geometries were calculated at the B3LYP/6-311G+(2d,p) level utilizing a “verytight” convergence criterion and an “ultrafine” integration grid. An anharmonic frequency calculation was performed to obtain the quartic and sextic centrifugal distortion constants for the normal isotopologue. Natural population analysis and Natural Bond Orbital (NBO)/Natural Resonance Theory (NRT) analysis of thiolium cations (C-atom protonated thiophenes) were performed with NBO 7.0.8.⁴⁵ 

Coupled-cluster calculations with single, double, and perturbative triple excitations [CCSD(T)] were performed using CFOUR.⁴⁶ Geometry optimizations were performed at the CCSD(T) level with the cc-pCVXZ basis set where X = D, T, Q, and 5. The structure computed at CCSD(T)/cc-pCVTZ was subsequently used for a second-order vibrational perturbation (VPT2) anharmonic frequency calculation by evaluating the cubic force constants using analytical second derivatives at displaced points.^47–49 The VPT2 calculation was used to obtain anharmonic frequencies for each isotopologue and provide vibration–rotation interaction constants (α_i), along with quartic and sextic distortion terms (in both the A and S reductions and I^r representation). The structure computed at the CCSD(T)/cc-pCVTZ level was also used for magnetic calculations to obtain the electron mass corrections to the rotational constants.⁵⁰ Additional computational corrections to the CCSD(T)/cc-pCV5Z structure are used to determine the best theoretical structure. These corrections were performed following the method prescribed by Heim et al.¹ and are specified below. The xrefit module of CFOUR was utilized to obtain the semi-experimental r_e^SE structure by least-squares fitting of the isotopologue equilibrium moments of inertia.

The synthetic approach to deuterium incorporation, via isotopic exchange, is ideal for rotational spectroscopy, enabling the broadband spectra of all isotopologues to be obtained in a single experiment. The characteristics of thiophene’s rotational spectrum cause the transitions of each isotopologue to be relatively distinct, thereby minimizing problems with the identification and assignment of transitions because of overlapping signals from multiple species. The ability to prepare and measure a single sample containing multiple species represents dramatic savings in time and effort, relative to an alternative approach involving the chemical synthesis and measurement of each isotopologue, individually.

Deuterium-containing isotopologues of thiophene were prepared in a manner analogous to the incorporation of deuterium in other five-membered heteroaromatic ring compounds and previous efforts on thiophene itself.^51–54 The acid-catalyzed reaction involved refluxing thiophene, D₂O, and fuming sulfuric acid at 100 °C for five days [Scheme 1(a)]. The duration of reflux may be varied to control the level of deuterium–hydrogen exchange. The five-day reaction time generated a mixture containing all nine of the possible deuterium-containing isotopologues of the normal isotopologue of thiophene. To shift the distribution of isotopologues toward those with a higher level of deuterium incorporation, the same reaction was conducted with refluxing extended over two weeks, with the addition of more H₂SO₄ and D₂O after seven days. This procedure yielded a small sample consisting of mostly tetra-deuterio thiophene, which provided transitions for the [2,3,4,5-²H, ³⁴S] isotopologue. While the extended time for reflux and replenishing the deuterium source did increase the amount of deuterium incorporation, this method is not sufficient for obtaining a sample that has only the tetra-deuterio isotopologue present.

The deuterium–hydrogen exchange reactions proceed via deuteriation/deprotonation at each carbon atom in the ring [Scheme 1(a)]. The preference for D/H exchange at C2 over C3 is consistent with the relative energies of the thiolium cation intermediates (Fig. 2). The 9.1 kcal/mol energy difference [B3LYP/6-311G+(2d,p)+ZPVE] is rationalized by delocalization of positive charge over more atoms via π conjugation when H/D⁺ is attached to C2, as shown in the accompanying resonance structures. A similar stabilization is not afforded when H/D⁺ is bonded to C3. This interpretation is drawn from the highest-contributing resonance structures from the natural resonance theory analysis and the natural population analysis charges of each thiolium ion, as displayed in Fig. 2. The differential rate of exchange leads to the dominance of [2-²H]-, [2,5-²H]-, and [2,3,5-²H]-thiophene in the rotational spectrum of the deuterium-enriched sample (after five days reflux), while the lowest-abundance thiophene-d_x species in the spectrum (excluding heavy-atom sulfur or carbon isotopologues) is [3,4-²H]-thiophene.

An alternate synthetic route was adapted to generate [3-²H]-thiophene⁵⁵ because that isotopologue cannot be prepared in high abundance using acid-catalyzed H/D exchange. A sample of [3-²H]-thiophene was generated from 3-bromothiophene by lithium–halogen exchange, followed by quenching with D₂O [Scheme 1(b)]. The resultant sample was primarily [3-²H]-thiophene, but the reaction did contain some [2-²H]-thiophene, along with other combinations of di- and tri-deuterio isotopologues (although in substantially smaller amounts than generated by acid-catalyzed H/D exchange). The rearranged and poly-deuteriated thiophenes presumably result from base-catalyzed H/D exchange in the strongly basic conditions. The reaction described in Scheme 1(b) allowed for the assignment of more transitions for the [3-²H] isotopologue and permitted the observation of sulfur and carbon heavy-atom isotopologues of [3-²H]-thiophene. Synthetic details and relevant NMR and mass spectra are included in the supplementary material.

The rotational spectrum of the ground vibrational state of the normal isotopologue of thiophene has been greatly expanded from previous works (5–40 GHz)^24,31–34 to include 40–360 GHz, along with some additional transitions between 12 and 18 GHz. The spectrum of thiophene, which includes the normal isotopologue and four heavy-atom isotopologues at natural abundance, is shown in Fig. 3. The distinct bands in the spectrum are due to ^aR_0,1 transitions that are degenerate with respect to K_c and decrease in J as they increase in K_a progressing toward higher frequency. As K_a continues to increase, these transitions lose degeneracy near K_a = 10 and eventually become degenerate with respect to K_a around K_a = 20. With the inclusion of transitions previously reported,^31,32,34 the dataset contains 3270 transitions with a J″ range from 1 to 93 and K_a″ range from 0 to 38; this is depicted in the data distribution plot provided in the supplementary material. Least-squares fitting of this large dataset requires a full sextic, distorted-rotor Hamiltonian (A and S reductions and I^r and III^r representations all employed) to obtain a low-error fit (σ = 0.022 MHz) with physically meaningful spectroscopic constants. The experimental κ value of −0.0917 indicates the highly asymmetric, barely prolate nature of thiophene and suggests the I^r representation to be somewhat more appropriate than III^r. The III^r representation, which is more appropriate to oblate asymmetric tops, has commonly been used in previous studies of thiophene. To allow for comparison with previous work, the A-reduced, I^r and III^r spectroscopic constants from this work are presented in Table I along with computational values and previously determined values from high-resolution IR studies.^18,20 The computed and experimental spectroscopic constants display excellent agreement. Compared to the experimentally determined spectroscopic constants (A-reduced, I^r representation), the computed rotational constants [CCSD(T)/cc-pCVTZ] are within 0.8%, the quartic distortion constants are within 2%, and the sextic constants are within 10%. This agreement provides excellent indication that the experimental constants (A reduction and I^r representation) are physically meaningful. Constants from the S-reduced Hamiltonian in the I^r representation are presented in Table II, and constants in the III^r representation least-squares fit are presented in Table SI. With the exception of H_K in the S reduction, I^r representation, all quartic and sextic constants were satisfactorily determined in both representations and in both reductions.

Clearly visible in the spectrum at natural abundance are two sulfur isotopologues ([³⁴S] and [³³S]) and two ¹³C isotopologues ([2-¹³C] and [3-¹³C]) of thiophene (Fig. 3). Transitions for each of these isotopologues from 40 to 360 GHz were combined with transitions from the literature (when within two times their reported experimental uncertainty)^25,34 and least-squares fit to a sextic, distorted-rotor Hamiltonian in the S reduction and I^r representation. For the ¹³C isotopologues (2.2% of the intensity of the normal isotopologue), over 960 distinct transitions were included in the least-squares fits. More than 1250 distinct transitions were least-squares fits for the [³⁴S] isotopologue (4.3% natural abundance), and over 250 transitions were least-squares fits for the much less abundant [³³S] isotopologue (0.75% natural abundance). No sextic distortion constants were able to be determined for the [³³S] isotopologue, unlike for the other three heavy-atom isotopologues, where each has some sextic distortion constants determined. The spectroscopic constants for all rare thiophene isotopologues are presented in Table II in the S reduction and I^r representation. All isotopologues have also been fit to a sextic, distorted-rotor Hamiltonian in the A reduction and I^r representation, which are presented in the supplementary material (Table SII), along with the data distribution plots for all other isotopologues (A reduction and I^r representation). Distortion constants that were not able to be experimentally determined were held constant to the predicted values [CCSD(T)/cc-pCVTZ]. The ³⁶S-isotopologue was sought using the rotational constants derived from the r_e^SE structure, but not identified due to its low natural abundance (0.01%).

Previously, rotational transitions of both mono-deuterio-, one di-deuterio ([3,4-²H])-, and tetra-deuteriothiophene isotopologues were reported.²⁴ The synthesis outlined in Scheme 1(a), after five days of reflux, resulted in a sample containing all possible mono-, di-, tri-, and tetra-deuterio isotopologues. A number of ³⁴S and ¹³C variants of these deuterium-enriched isotopologues were also observed at their natural abundances. In total, 15 isotopologues of thiophene, including 11 that were previously unreported, were observed from this sample. For each isotopologue, rotational transitions were least-squares fits to a sextic, distorted-rotor Hamiltonian in the A and S reductions and the I^r representation (S reduction in Table II and A reduction in the supplementary material, Table SII). The deuterio isotopologues follow a similar spectral pattern to the normal isotopologue, with strong ^aR_0,1 bands; a portion of the spectrum of the deuterium-containing isotopologues is displayed in Fig. 4. The band structure of thiophene is such that the spacing between bands and the transitions within the bands allows for clear visibility of each of the various isotopologues. While the [2,3,4,5-²H, ³⁴S] species was observed in the original five-day deuterium-enriched sample, the measured transitions reported in Table II come from the sample containing primarily tetra-deuterio thiophene (14-day sample). Transitions from previous microwave work²⁴ of the [2-²H] and [3-²H] isotopologues that agreed to within 2σ of their experimental uncertainty were included in their respective datasets (data distribution plots provided in the supplementary material).

In the [3-²H]-enriched sample [Scheme 1(b)], [3-²H, ³⁴S]-thiophene was the most abundant heavy-atom isotopologue, and its enhanced intensity allowed for the measurement of 664 distinct transitions. The presence of deuterium at the three-position breaks the C_2v symmetry, which results in four unique mono-¹³C isotopologues (1.1% the intensity of the [3-²H] isotopologue). Over 100 transitions were measured and least-squares fit for each [3-²H, ¹³C] isotopologue. The incomplete specificity of the reaction also resulted in deuterium incorporation at the two and five positions and provided additional transitions for the [2-²H, ³⁴S] isotopologue. Inasmuch as ³³S is the least naturally abundant isotope observed, only 57 transitions for [3-²H, ³³S]-thiophene were able to be included in its least-squares fit. Like with the [2,3,4,5-²H, ³⁴S] isotopologue, the small number of transitions and relatively large error in A₀ and B₀ means that the spectroscopic constants of [3-²H, ³³S] are not confidently determined. Therefore, neither of these isotopologues were included in the determination of the structure.

The semi-experimental equilibrium structural parameters (r_e^SE) of thiophene were determined by least-squares fitting of the equilibrium moments of inertia. The experimental rotational constants (B₀) of 24 isotopologues provide 72 moments of inertia to determine eight independent structural parameters (due to the C_2v symmetry of thiophene). The B₀ constants from both the A- and S-reductions (I^r representations) were converted to their determinable constants (B″) to remove effects of centrifugal distortion using Eqs. (S1)–(S6) (see the supplementary material)⁵⁶ and then averaged. The very close agreement of the B″ constants (Table SIII) confirms that each isotopologue transition dataset is well treated by both the A- and S-reduced Hamiltonians, which produce physically meaningful spectroscopic constants. The averaged B″ constants from the A and S reductions were further corrected to their equilibrium rotational constants by applying the computational vibration–rotation corrections (half the sum of the vibration–rotation interaction constants, α_i) and the computational electron mass corrections, as described in Eq. (1).^1,30 The electron mass corrections are obtained from the diagonal elements of the electronic contribution to the rotational g-tensor, g^ββ, multiplied by the electron–proton mass ratio (m_e/M_P) and the computed equilibrium rotational constant [B^β_CCSD(T)], where β = a, b, or c,^1,30

The deviation in the experimentally determined inertial defect (Δ_i = I_c − I_a − I_b) from the ideal inertial defect value of zero for a planar molecule is much reduced by the inclusion of these computational corrections. The effect of each of the corrections on the inertial defect is shown in Table III for each isotopologue included in the structure. Table III includes Δ_i r₀, the inertial defects from only the experimental B″ values relating to an r₀ structure, Δ_i r_e^SE (vib. corr. only), the inertial defects from the B″ values applying the vibration–rotation corrections only, and Δ_i r_e^SE, the inertial defects calculated from the B″ values applying the vibration–rotation and electron mass corrections. The averages (⁠$x̄$⁠) of these respective inertial defects, as well as the individual values for each isotopologue, indicate that both the vibrational and electronic corrections are needed to bring Δ_i closer to the ideal value of zero and obtain the best possible equilibrium constants (B_e). It is further clear that, including both corrections, the inertial defects are highly systematic and fairly constant. The low sample standard deviations (s) strongly indicate the high level of precision in the equilibrium moments of inertia used to calculate Δ_i r_e^SE and that the remaining deviations from zero are due to a yet unidentified systematic effect. Any random errors in the rotational constants (including in the vibration–rotation corrections) larger than a few tens of kHz would cause the value of s to be greater than what is observed. The ratio of $x̄$/s in the last column (Δ_i with both corrections applied) is ∼7. It is also notable that all the shifts due to the electron mass correction are +0.010 30 or +0.010 31, even though the corresponding shifts in the individual rotational constants vary with isotopologue. The constancy in this shift is, of course, reflected in the standard deviation being the same in the two corresponding columns.

The averaged B″ values and computational corrections for each isotopologue are submitted to the xrefit module in CFOUR,⁴⁶ which least-squares fits the semi-experimental equilibrium moments of inertia and determines the structural parameters. The r_e^SE structural parameters are presented in Fig. 5(a) and in Table IV, alongside the computed structural parameters [CCSD(T)/cc-pCV5Z] and experimental structural parameters determined previously.^25,27,28,30 Figure 5(b) presents the number of isotopologues with substitution(s) at the designated position used in the least-squares fit of the structure from 24 isotopologues. (For example, 3 at the C2 position indicates that its substitution is present in three isotopologues: [2-¹³C], [2,5-²H, 2-¹³C], and [3-²H, 2-¹³C].) The structural parameters are presented graphically in Fig. 6.

To provide the best comparison to the semi-experimental equilibrium structure, the following four computational corrections to the CCSD(T)/cc-pCV5Z structure have been implemented to address residual error in the computational structure:

1. Residual basis set effects beyond cc-pCV5Z.

2. Residual electron correlation effects beyond the CCSD(T) treatment.

3. Effects of scalar (mass–velocity and Darwin) relativistic effects.

4. The diagonal Born–Oppenheimer correction (DBOC).

These four corrections are obtained in a manner similar to pyrimidine¹ using Eqs. (2)–(8), as described below.

1. To estimate the correction needed to approach the infinite basis set limit, equilibrium structural parameters obtained with the cc-pCVXZ (X = T, Q, and 5) basis sets were extrapolated using the empirical exponential⁵⁸ 

   $Rx=R∞+Ae−Bx,$

   (2)

   where R(x) are the values of the parameters obtained using the largest three basis sets (x = 3, 4, and 5) and R(∞) is the desired basis set limit estimate. Using these three basis sets, Eq. (2) can be expressed as follows:

   $R∞=−R42−R3R5R3+R5−2R4.$

   (3)

   The correction to the structure due to a finite basis set is then estimated by the following equation:

   $ΔRbasis=R∞−RCCSD(T)/cc-pCV5Z.$

   (4)
2. Residual correlation effects are assessed by optimizing geometries at the CCSDT(Q) level⁵⁹ and then estimating the correlation correction as in the following equation:

   $ΔRcor=RCCSDT(Q)−RCCSD(T).$

   (5)

   Due to the expense of calculations at the CCSDT(Q) level, the cc-pVDZ basis in the frozen-core approximation is used for these calculations.
3. Relativistic corrections are obtained by subtraction of the parameters with a standard non-relativistic calculation from the equilibrium structure obtained with the X2C-1e variant of coupled-cluster theory,^60–62 as in the following equation:

   $ΔRrel=RCCSD(T)/cc-pCVTZSFX2C−1e−RCCSD(T)/cc-pCVTZNR.$

   (6)
4. The diagonal Born–Oppenheimer correction (DBOC)^63,64 is obtained from the following equation:

   $ΔRDBOC=RSCF/cc-pVTZDBOC−RSCF/cc-pVTZNR.$

   (7)

   Here, the first value is obtained by minimizing the DBOC-corrected SCF energy with respect to nuclear positions, and the latter is again the traditional calculation.

The sum of these corrections is used to obtain the best equilibrium structural parameters, given by the following equation:

The values of each correction, ΔR(best), are applied to the CCSD(T)/cc-pCV5Z structural parameters and are presented in Table V and Fig. 6. We will refer to them as the “best theoretical estimate” (BTE) parameters.

In the r_e^SE structure of thiophene, the largest computed statistical (2σ) bond distance uncertainty is in R_C2–C3 (0.000 31 Å), and the largest bond angle uncertainty is in A_S1–C2–H (0.028°). Exactly how reliable these error estimates are, however, is a more subtle question. Comparing current r_e^SE to the previously published r_e^SE structure (Table IV), which applied vibration–rotation and electron mass corrections computed using DFT methods³⁰ to the rotational constants of eight isotopologues,²⁵ none of the previous parameters fall within the 2σ error bars of our current r_e^SE structure except for the R_C2–H distance. For most parameters, the earlier r_e^SE structure³⁰ is not in close agreement with our BTE structure, either. The differences between the two r_e^SE structures are presumably due to a combination of the incorporation of the precisely determined moments of inertia for more isotopologues and to differences in computational treatment.

The large number of isotopologues incorporated in our current studies poses some new challenges with respect to analysis and interpretation of data. Therefore, our group developed a new method to analyze the contribution and importance of each isotopologue to the r_e^SE structure determination. A traditional substitution structure (r_s) from a Kraitchman analysis²⁶ relies on single isotopic substitution, which is sufficient to determine the structure.²⁵ Originally, this calculation employed ground-state constants (B₀), but it can also be applied to equilibrium rotational constants. The success of the Kraitchman analysis demonstrates that there is, at least, in principle, sufficient information content in the value of the rotational constants for the normal form and the singly substituted forms to locate every atom in the molecule. We have found, however, that the use of multiply substituted isotopologues provides useful additional structural information.^1,38 The parameters become highly over-determined, and the experimental uncertainty can be pushed well below the value achievable with smaller datasets, i.e., single substitutions. The analysis of isotopologue incorporation was performed with an iterative analysis, where we began with the normal isotopologue and five mono-substituted isotopologues ([2-¹³C], [3-¹³C], [³⁴S], [2-²H], and [3-²H]) as the “core set” and then added each of the remaining isotopologues, one-at-a-time, to the least-squares fit of the structure in xrefit. The value of the square root of the sum of the squared relative errors, which we call the δ (r_e^SE), was determined for each isotopologue addition by using the following equation:

where P is each structural parameter determined by xrefit and σ_p is the respective (1σ) error of each parameter. The isotopologue that results in the lowest value of the δ (r_e^SE) for the extended set is then incorporated into the structure determination for a total of seven isotopologues. The procedure repeats again, testing each remaining isotopologue and incorporating the one that gives the smallest δ (r_e^SE) into the fit, for a new total of eight isotopologues. This procedure repeats until all 24 isotopologues are incorporated. As this iterative process could be tedious, a program, xrefiteration, has been written to automate it. The code will be discussed in more detail in a publication on the improved semi-experimental equilibrium structure determination of pyridazine.⁶⁵ The resultant δ(r_e^SE) values from each iteration are plotted in Fig. 7, along with δ(r_e^SE) for the bond distances and the angles individually.

Upon examination of the xrefiteration results, we observe that the addition of the [2,3-²H] and [2,4-²H] isotopologues (the last two) causes δ(r_e^SE) with 24 isotopologues to increase to nearly the same as δ(r_e^SE) with 12 isotopologues. In an attempt to understand this behavior, we excluded data from the last two isotopologues ([2,3-²H] and [2,4-²H]) and determined a new semi-experimental equilibrium structure, labeled r_e^SE (r₂₂), which is provided in the supplementary material (Table SIV) alongside the r_e^SE structure. Despite the modest increase in the uncertainties obtained by including the last two isotopologues, the structural parameters themselves [r_e^SE vs r_e^SE (r₂₂)] are hardly changed at all. The most substantial impact of the final two isotopologues, beyond raising δ(r_e^SE), is to move the r_e^SE values of R_C2–H and A_S1–C2–H closer to their BTE values (Fig. 8). This suggests that these isotopologues are providing important structural information about the locations of C2 (and C5) that is not contained in the isotopologue dataset from the previous 22 isotopologues.

A comparison of the various correction terms for the best theoretical estimate of the equilibrium structure for thiophene (Table V) to equivalent data for pyrimidine¹ shows that values for thiophene are somewhat larger, but not by a great deal. This leads us to suspect that the BTE results for thiophene might be of comparable accuracy to those obtained previously for pyrimidine. Referring to the last two columns of Table V or to Figs. 6 and 8, we observe that, for several structural parameters of thiophene, the BTE values are within 2σ of the r_e^SE values, or nearly so, and the corresponding residuals are similar in magnitude to those found in pyrimidine. In contrast to pyrimidine, Fig. 8 shows that for two of the bond distances (R_C2–C3 and R_C2–H) and one of the angles in thiophene (A_S1–C2–H), the BTE values are well outside 2σ of their corresponding r_e^SE values. It is easily seen that all the offending parameters involve C2 and the H atom on that carbon. We suspect that this problem arises because these two atoms lie close to the b principal axis, a classically problematic issue in using moments of inertia to determine structures.²⁵ In a planar molecule, the atom a axis coordinates determine I_b and the b axis coordinates determine I_a. An atom being near the b axis will hinder determination of its a axis coordinate and conversely. The uncertainties in atom coordinates are typically magnified when those atoms lie near any principal axis, but not on it. Sometimes, the principal axes can be substantially rotated by isotopic substitution to move such atoms away from the axis and greatly mitigate this problem. Not much rotation can occur in thiophene, however, because the heavy sulfur atom must always be very near the a axis.

To obtain a better understanding of the origin of the structural errors in thiophene, we expressed the various structure types in a new coordinate system shown in Table VI. We use rectangular coordinates aligned with the symmetry axis, but with the sulfur atom at the origin. There are four unique a axis coordinates and four unique b axis coordinates. We find that the agreement between BTE and r_e^SE for all b axis coordinates is excellent. The a axis coordinate residuals for C3 and H3 are fairly small (especially considering that the values themselves are large), but the a axis coordinate residuals for C2 (0.000 52 Å) and, especially, H2 (−0.002 63 Å) are clearly where the problem lies. Because there is no reason to believe that the BTE value of the a axis coordinates would be less accurate than that of the b axis coordinates, it seems reasonable that any accuracy problems in the present thiophene r_e^SE determination trace back to the proximity of the C2 and H2 atoms to the b principal axis, for which we have not yet found a satisfactory solution.

The current study explores the limits of how accurately and precisely an r_e^SE structure can be determined by rotational spectroscopy with CCSD(T) corrections. Despite obtaining an r_e^SE structure for thiophene that is a substantial improvement on previous structure determinations for molecules containing sulfur, the agreement found for some parameters is somewhat problematic compared to the outstanding agreement found for r_e^SE of pyrimidine. For pyrimidine, all structural parameters of best theoretical CCSD(T) r_e are within the statistical uncertainties of semi-experimental r_e^SE with CCSD(T) corrections. For thiophene, this is not the case, despite a substantial over-determination of its geometry with moments of inertia from 24 isotopologues. Our analysis indicates that the proximity of C2 and H2 to the b principal axis is the source of the problem, preventing the locations of those atoms from being determined to sufficient accuracy. It was not possible in this study to overcome the problem by the inclusion of isotopologues that cause a large rotation of the principal axes. The large mass of the sulfur atom in thiophene prevents any isotopic substitution from causing the necessary large rotation of the principal axes. Given the accuracy of the determination of all atom positions except those near the b axis, this study provides strong evidence that satisfactory agreement between r_e^SE and the best computed parameters is as likely for a molecule containing a third-row atom, as for the one which does not. Future investigations of other molecules with third-row atoms (Si, P, S, Cl, etc.), but without any atoms lying near a principal axis, would provide an opportunity to confirm this assertion.

See the supplementary material for least-squares fitting files for all thiophene isotopologues, spectroscopic constants in the A-reduction, data distribution plots for isotopologues, determinable constant equations, determinable constants from the A and S reductions, output files from computations, output files from the structure fitting programs (xrefit and xrefiteration), synthetic details, and characterization data.

We acknowledge funding from the National Science Foundation (Grant Nos. RJM CHE-1664912, CHE-1954270, and JFS CHE-1664325). We also acknowledge funding from the National Institutes of Health for the support of the shared UW-Madison departmental mass spectrometry instruments (Grant No. NIH 1S10 OD020022-1). K.K. thanks KAKENHI, Grants in-Aid for Scientific Research by the Ministry of Education, Culture, Sports, Science, and Technology of Japan (Grant Nos. 26108507, 15H03646, and 19H05074), and the National Astronomical Observatory of Japan. We acknowledge Michael McCarthy for the loan of an amplification multiplication chain. We thank the Harvey Spangler Award (to B.J.E.) for the funding that supported the purchase of a zero-bias detector. We thank Maria A. Zdanovskaia and Andrew N. Owen for their assistance and thoughtful conversations.

The data that support the findings of this study are available within the article and its supplementary material.