Millimeter-Wave Spectroscopy of Methylfuran Isomers: Local vs. Global Treatments of the Internal Rotation

Methylfurans are methylated aromatic heterocyclic volatile organic compounds and primary or secondary pollutants in the atmosphere due to their capability to form secondary organic aerosols in presence of atmospheric oxidants. There is therefore a significant interest to monitor these molecules in the gas phase. High resolution spectroscopic studies of methylated furan compounds are generally limited to pure rotational spectroscopy in the vibrational ground state. This lack of results might be explained by the difficulties arisen from the internal rotation of the methyl group inducing non-trivial patterns in the rotational spectra. In this study, we discuss the benefits to assign the mm-wave rotational-torsional spectra of methylfuran with the global approach of the BELGI-Cs code compared to local approaches such as XIAM and ERHAM. The global approach reproduces the observed rotational lines of 2-methylfuran and 3-methylfuran in the mm-wave region at the experimental accuracy for the ground vt=0 and the first torsional vt=1 states with a unique set of molecular parameters. In addition, the V3 and V6 parameters describing the internal rotation potential barrier may be determined with a high degree of accuracy with the global approach. Finally, a discussion with other heterocyclic compounds enables the study of the influence of the electronic environment on the hindered rotation of the methyl group.


Introduction
Furan and its derivatives (furans) are heterocyclic organic compounds belonging to the family of oxygenated five-membered aromatic molecules. They are primary and secondary pollutants in the atmosphere emitted from multiple sources such as oil refining, coal mining and gasification, biomass and fossil fuel or waste combustions [1,2]. Furans and other heterocyclic compounds such as pyrroles (with N as the heteroatom in the ring) are produced by the pyrolysis of cellulose and are major components of the emission of wildfire burnings [3]. Such compounds as the volatile furan (C 4 H 4 O) and methylfurans (C 5 H 6 O) (see Figure 1) or the less volatile furaldehyde (C 5 H 4 O 2 ) have emission levels 70 to 120 times higher compared to CO [4]. Once furans are emitted, they will undergo gas phase chemistry and, to an extent, will be photolyzed at actinic wavelengths to produce tropospheric ozone or they will react with the main atmospheric oxidants (OH, Cl atoms, ozone) leading to the formation of secondary organic aerosols (SOA) [5]. A similar process occurs during the night by the oxidation of furans with NO 3 radicals [1][2][3]. These SOA affect the climate both directly and indirectly and large uncertainties still remain on their radiative forcing [6]. For all these reasons, there is a fundamental interest to monitor furans directly in the atmosphere or in atmospheric simulation chambers to characterize their reactivity and their ability to produce SOA. Up today, the monitoring of volatile organic compounds (VOC) involved in biomass burning emissions such as furans is mostly performed off-line with chemical methods based on gas chromatography (GC) and mass spectrometry (MS). A GC-MS instrument provides extensive chemical details of discrete air samples collected in wildfires or during a controlled combustion process in the laboratory. These measurements can be complemented by real-time middle or low resolution measurements with an open-path Fourier transform infrared spectroscopy instrument [4,7]. To the best of our knowledge, high-resolution spectroscopic methods were never used to measure on-line individual signatures of furans at trace level in realistic environment. Trace gas detection in complex chemical mixtures with IR and THz techniques are generally focused on a limited panel of 55 light atmospheric molecules listed in the high-resolution HITRAN database [8]. No high resolution data of furans are present in this database, only vibrational cross-sections of furan, methylfurans (MF) and furaldehyde may be used for quantitative spectroscopy in the IR at middle/low resolution. This might be explained by the difficulty to resolve highly congested room temperature rovibrational spectra of these molecules. To overcome this challenge, one first need to measure the spectra in laboratory and assign pure rotational and low-frequency rovibrational spectra of furans at room temperature. Several groups have performed high-resolution pure rotational spectroscopic studies of furans in millimeter-wave at room temperature. We can cite the furan ground state analysis of Wlodarczack et al. [9]; the ground and first torsional states of trans-furfural and ground and excited state of cisfurfural analyses of Motiyenko et al. [10] assigning in total seventeen states. In 1971 Ogata and Kozima recorded the microwave spectrum of the 3-MF between 8.2 and 30 GHz and assigned 20 A and 16 E species transitions [11] and more recently, the ground and first torsional states of the 2-MF isomer were analyzed by Finneran et al. requiring a specific treatment of the internal rotation of the methyl group [12]. Finally, we can mention the study of furan performed by Tokaryk et al. using FT-Far-IR spectroscopy based on a synchrotron source which constitutes, to the best of our knowledge, the only high-resolution rovibrational study of a furan compound [13].
In this study, new high-resolution millimeter-wave measurements have been performed on 3-MF with a versatile solid-state spectrometer already used for the investigation of other VOC such as methoxyphenol [14] or catechol [15] and gas phase monitoring in an atmospheric simulation chamber [16,17]. As in Ref. [12] with 2-MF, the spectroscopic assignment was complicated by the large amplitude motion associated with the internal rotation of the methyl group. Unlike the work of Finneran et al. where the internal rotation was treated with two local approaches, we have succeeded in this study, using a global treatment, to reproduce the spectra in the ground and first excited torsional states at almost the experimental accuracy with a common set of molecular parameters. Amongst them, the V 3 and V 6 constants of 2-MF and 3-MF are obtained with a high degree of accuracy allowing to characterise finely the internal rotation potential for both MF isomers.

Results
The optimized geometries of 2-MF and 3-MF isomers at the B3LYP/aug-cc-pVTZ level of theory are presented in Figure 1. Both isomers are prolate asymmetric tops with Ray's asymmetry parameter κ 2−MF = −0.69 and κ 3−MF = −0.73. Their permanent dipole moment is oriented in the (a, b) plane with the |µ a |/|µ b | ratio of around 0.6 [18] and 1.9 [11], respectively. Both a-type and b-type transitions are observed in the MF rotational spectra: the most intense transitions are expected to be b-type for 2-MF and a-type for 3-MF. One difficulty in the assignment process, already mentioned by Finneran et al. in Ref. [12], lies in the possible mixing of energy levels arising from comparable asymmetry and internal rotor splitting. Therefore, we were careful to keep only consistent labels for the energy levels. This difficulty which is generally observed in the rotational-torsional spectra of other small organic molecules containing a methyl internal rotor will be discussed more in detail in Section 3.1.

Global Fit of the Vibrational Ground State and First Torsional State of 2-MF
Finneran et al. [12] recorded in 2012 a broad spectrum of 2-MF from 8.7 GHz to 960 GHz. In the 8.7-18.3 GHz range they used a room-temperature chirped-pulse Fourier transform microwave spectrometer. Above 75 GHz, the mm-wave and submm-wave rotational spectra of 2-MF have also been measured by those authors [12] using a direct absorption flow cell spectrometer [19]. The spectra have been obtained in the following spectral ranges: 75-125 GHz, 130-450 GHz, 460-650 GHz, 700-732 GHz, and 760-969 GHz. In the present work, we used a similar spectrometer with an updated technology, a doublepass absorption cell and a room temperature detection based on zero-biased detector (ZBD) [20]. As a test, we decided therefore to measure again the 2-MF spectrum with the first stage of multiplication (×9) covering the 70-110 GHz range. In Figure 2, a part of this spectrum is presented from 86.4 GHz to 86.9 GHz, and we tagged the rotational lines previously assigned in Ref. [12].
Only a few lines remain unassigned. As already mentioned by Finneran et al. [12], apart from the methyl torsion band at 158 cm −1 , two other low frequency modes occur below 400 cm −1 associated with an out-of-plane bending (237 cm −1 ) and an in-plane bending (339 cm −1 ) of the methyl group. Those states are most likely populated at room temperature and thus can be responsible for those unassigned lines. In Ref. [12], Finneran et al. have used the eXtended Internal Axis Method (XIAM) [21] and the Effective Rotationaltorsional Hamiltonian (ERHAM) [22] to fit separately the ground state (v t = 0) and the first excited torsional state (v t = 1) parameters within two "local" approaches. Unlike the XIAM code which failed to fit v t = 1, the ERHAM code succeeded in fitting both v t = 0 and v t = 1 states separately within their experimental accuracy (see Table 1). We plotted in Figure 2b the simulated spectrum with the ERHAM parameters listed in Table 1 of [12] using the prediction published in their supplementary data. Most of the observed rotational lines are predicted taking into account assignments from v t = 0 and v t = 1 states. The remaining unassigned lines belonging to other excited states involving either v t > 1 states or higher frequency vibrational modes are significantly less intense. In this work, we decided to perform a new fit using a global approach with a unique set of parameters for both the v t = 0 and v t = 1 states by means of the BELGI-C s code [23]. GHz frequency range. Blue and red sticks correspond respectively to assigned v t = 0 and v t = 1 lines in Ref. [12]. (b) Simulated spectrum with the v t = 0 (blue) and v t = 1 (red) parameters fitted in Ref. [12] using ERHAM (local approach). (c) Simulated spectrum of the present work with the v t = 0 (blue) and v t = 1 (red) using BELGI-C s (global approach) and the fitted parameters from Table S1. For XIAM, as shown in Ref. [12], due to the lack of sufficient higher order distortion constants, as well as high order torsion-rotation terms, their fit leads to several difficulties at higher K a values. For this reason, we decided to start the global fit with BELGI-C s using only the assigned transitions from ERHAM for v t = 0 and v t = 1 states. However, the conversion of the quantum numbers used to label energy levels and molecular transitions in ERHAM to the labelling scheme used in the BELGI-C s is not straightforward since the values of K a and K c are approximately determined by the (assumed) order of energy levels for rigid asymmetric rotors. These labels may be incorrect when serious rotational-torsional interactions occur. Thus, for the global fit with BELGI-C s many of the K c values were reassigned. In the range up to J = 70 (limit of the BELGI-C s code), only 14 transitions (which correspond to 7 distinct frequencies) with K a = 53 and 54 used in ERHAM were excluded in our BELGI-C s fit, as we could not find correct labelling for those levels. As for the ERHAM fit, the BELGI-C s fit succeeded to reproduce the experimental 2-MF linelist of Ref. [12] at almost their experimental accuracy.
The complete set of the Rho Axis Method (RAM) parameters including centrifugal distortion and higher order terms used in the BELGI-C s is provided in the supplementary materials (Table S1). In Table 1, we present the fitted rotational constants for both the (v t = 0) and (v t = 1) states, after transforming them into the principal axis method (PAM). In addition, the global approach allows us to fit the height of the internal rotation barrier V 3 and its higher-order correction V 6 with a high degree of accuracy. This is also the case for some structural parameters like the angles between the internal rotation axis and the principal inertia axes. We notice the good agreement of these experimentally derived parameters with those obtained from the B3LYP/aug-cc-pVTZ calculations presented in Section 4.1.1 and discussed in Section 3.2. Table 1. Torsion-rotation constants of 2-MF for the v t = 0 and v t = 1 states obtained from DFT calculations (this work), from the fits of microwave and mm-wave data with local XIAM and ERHAM approaches (Ref. [12]) and global BELGI-C s approach (this work). The complete set of parameters obtained for XIAM and ERHAM can be found in Tables 1 and 2 of Ref. [12]. The complete set of parameters obtained for BELGI-C s is given in the supplementary Table S1. The numbers in parentheses represent 1 σ standard deviation in units of the last significant digit. All parameters are in the Principal Axis System (PAS).
Ref. [12] Ref. [12] Ref. [  a Calculation with the aug-cc-pVTZ basis set. b The complete set of higher order tunneling terms used in ERHAM are listed in Tables 1 and 2 of Ref. [12]. c Parameters in the PAS after transforming them from the RAM. The complete set of RAM parameters including centrifugal distortion and higher order terms used in our fit are listed in Table S1 of the supplementary data. d In XIAM F is a derived parameter from F 0 , see supplementary data in Ref. [12]. In BELGI-C s F was fixed to the best value obtained by fitting it. e Standard error of the analytical fit based on Equation (5) as given by Gnuplot [24]. f Reduced barrier height: s = 4V 3 /(9F). g N corresponds to the number of assigned transitions for v t = 0 and v t = 1. In the BELGI-C s fit 14152 and 3742 correspond to the number of transitions fitted. h In XIAM only 8006 unblended lines were included in the fit. i In the ERHAM fit 11,793 and 2580 correspond to the number of distinct frequencies used for the fits of Ref. [12]. j WRMS is the weighted unitless root mean square deviation of the fit: 1 ) and calculated ( f calc. i ) frequencies and ∆ i the estimated error on the frequency measurement. k RMS is the Root mean square deviation of the fit in kHz. Figure 2c shows 500 MHz of a 2-MF spectrum simulated with the unique set of molecular constants obtained from the BELGI-C s fit (The complete line list is provided in Table S3 of the supplementary material). As for Figure 2b, a simulated spectrum based on the ERHAM fit, the rotational lines belonging to the v t = 0 and v t = 1 states are plotted respectively in blue and in red. A direct comparison between the local ERHAM and global BELGI-C s fits may be performed. Both predictions are able to reproduce reasonably the measured line positions of the 2-MF spectrum. Nevertheless, some differences are observed between these two simulated spectra especially when considering the relative intensities. Indeed the Figure 2b,c exhibit a few line intensities which are not in agreement with the experimental ones (e.g: the v t = 1 line centered around 86,620 MHz for the ERHAM prediction or the v t = 0 line centered around 86,786 MHz for the BELGI-C s prediction). The differences between local (XIAM and ERHAM) and global (BELGI-C s ) are discussed in more details in Section 4.1.2.

Global Fit of the Vibrational Ground State and First Torsional State of 3-MF
Concerning 3-MF, a microwave study was performed in 1971 by Ogata et al. [11] in the 8.2-30 GHz range, providing a first estimation of the rotational constants for the vibrational ground state. This study allowed them to determine the permanent dipole moment from Stark measurements (1.03 D) and constants associated to the CH 3 internal rotation like the V 3 barrier height (381 cm −1 ) from the observed A-E splittings of 16 transitions. In our present work the mm-wave spectral analysis of 3-MF allowed us to identify 2463 new transitions in the v t = 0 state and 2017 in the first excited v t = 1 state up to J = 50.
As shown in Figure 3, a simulated spectrum (red) consisting of the v t = 0 A-(green) and E-species (blue) as well as v t = 1 A-(magenta) and E-species (brown) rotational lines predicted with our BELGI-C s fit reproduces the experimental spectrum (black) well. For the frequency range of 70 -220 GHz, J max = 50 and K a,max = 31 for v t = 0 and v t = 1 was sufficient to assign the most intense lines. Remaining unassigned weaker lines are probably associated to the out-of-plane and in-plane methyl bending modes predicted respectively at 267 cm −1 and 321 cm −1 with anharmonic B3LYP/aug-cc-pVTZ frequency calculation. To compare the efficiency of local and global approaches (see Section 3.1), the same data set (excluding one line) was fitted using the XIAM code for the v t = 0 and v t = 1 states separately. The results of these fits are given in Table 2 along with the theoretical B3LYP/aug-cc-pVTZ calculation. Again, as for the 2-MF, our PAM values fitted by the BELGI-C s code agree reasonably well with the theoretical calculations. The complete set of constants (19 parameters) obtained from the BELGI-C s global fit in the RAM axis system is presented in Table S2 and the line list deduced from our global BELGI-C s fit is given in Table S4 of the supplementary materials. Table 2. Torsion-rotation constants of the 3-MF obtained from the mm-wave (70-220 GHz) analysis using the XIAM code for v t = 0 and v t = 1 separately, and BELGI-C s with a global fit analysis. These values are compared with quantum chemistry calculations (B3LYP method) and with the microwave analysis performed by Ogata et al. [11]. The numbers in parentheses represent 1 σ standard deviation in units of the last significant digit. All parameters are in the principal axes system.  (5) as given by Gnuplot [24]. g Reduced barrier height: s = 4V 3 /(9 × F) h N corresponds to the number of assigned transitions for v t = 0 and v t = 1. i WRMS is the weighted unitless root mean square deviation of the fit: ) and calculated ( f calc. i ) frequencies, for v t = 0 N b = 1845, for v t = 1 N b = 1538 and ∆ i is the measurement uncertainty. j RMS is the root mean square deviation of the fit in kHz. k Previous fit of Ogata et al. [11] using only A species transitions and refined from A-E splitting analysis with A A − A = 2.33 MHz and B A − B = 0.00 MHz.

"Local" versus "Global" Approach
A "local" method treats separately each torsional v t state, whereas a "global" method explicitly takes into account all the interactions between the different v t states. Both XIAM and ERHAM programs adjust the parameters of each torsional state separately: they belong to "local" approaches. BELGI-C s on the other hand allows, for a given vibrational state, the fit of several torsional states together within the same rotation-torsion Hamiltonian matrix, see Section 4.1.2. A comparison of the results obtained with a "local" approach with those issued from a "global" approach has already been introduced in several papers [23,25,26]. Local approaches have shown their efficiency and high speed, especially in the first steps of assignments, but they sometimes fail when the potential barrier is low or when reaching higher excited torsional states. Finally one of the differences between XIAM and ERHAM is that the potential barrier does not appear explicitly in the ERHAM method. The information about the barrier is thus hidden in the tunneling coefficients of the internal rotation energy.
For 2-MF, as mentioned by Finneran et al. [12], ERHAM computes very quickly and accurately the line positions and intensities up to J = 120 and provides much better RMS deviations (108 kHz for v t = 0 and 113 kHz for v t = 1) than XIAM (263 kHz and 16 MHz respectively). ERHAM used five and ten tunneling parameters, (the Bkpr q parameters listed in Tables 1 and 2 of Ref. [12] and defined in Ref. [22,27]) to fit v t = 0 and 1 respectively. The XIAM fit stays of rather poor quality, even after they added four sextic centrifugal distortion, and three rotation-torsion parameters, (see third column of Table 1 [12]). On the other hand, XIAM provides a direct value for the potential barrier [21].
Our global fit for 2-MF with BELGI-C s simultaneously includes both the v t = 0 and 1 rotational transitions and achieves a similar quality as ERHAM. The reserved memory by the BELGI-C s program allows the construction of the Hamiltonian matrix for J up to 70 with an Hamiltonian matrix dimension is (2J + 1) × 9. Even though we only assigned transitions of v t = 0 and 1, we take into account interactions within the lowest nine torsional states, a truncation usually allowing for the best lowest states energy calculation. The RMS deviations obtained by BELGI-C s are 107 and 114 kHz for v t = 0 and 1 respectively. As shown in our supplementary Table S1, we use in total 21 parameters (plus F which is fixed): the three rotational constants A, B and C, the off diagonal inertia moment, D ab , the five fourth order centrifugal distortion terms, and only one sextic centrifugal constant (H K ). We also fit V 3 and V 6 , ρ the coupling term between internal and global rotation, seven fourth order and one six order rotation-torsion terms.
Compared to the XIAM values, our values, after transforming them into the PAM, are in reasonable agreement, considering the difference between the two methods (see Table 1). One of the differences between XIAM and BELGI-C s is that XIAM fits the three angles ∠(i, x), ∠(i, y) and ∠(i, z) between the internal axis of the methyl group and the principal axes, whereas BELGI-C s fits the D ab term. We have to diagonalize the inertia tensor, using a procedure described in [25], to get the eigenvalues which are the A, B and C values in the PAM. The eigenvectors generated by this diagonalization of the tensor of inertia are related to the angles ∠(i, x), ∠(i, y) and ∠(i, z) reported in Table 1. Compared to XIAM our A, B, C parameters in the PAM only differ by 0.0002, 0.05, 0.01 % to XIAM. The ∠(i, x), ∠(i, y) and ∠(i, z) angles have also close values, showing that XIAM does well with determining the structure, even with a poor quality of fit. Compared to the theoretical B3LYP calculation, the agreement of our PAM values A, B, C from the BELGI-C s code is 0.9, 0.2, 0.4 %. The V 3 from BELGI-C s only differs by 1.8 % with its value from XIAM, and by 11.8 % with the B3LYP calculation, which is also satisfying. We note that a good order of magnitude is obtained as well for V 6 by the theoretical method compared to the BELGI-C s value. The angle between the PAM axis and the RAM z axis is listed in Table 1 as "β" to be coherent with ERHAM notation. We called this angle θ RAM in our previous paper [28] (see Equation (6) of the theoretical section). Its order of magnitude is conserved between BELGI-C s and ERHAM, but ERHAM has two distinct values for v t = 0 and 1. Finally we note that in Table 1, we do not indicate centrifugal distortion terms for BELGI-C s . Indeed, their values are fitted in the RAM system, but there is no way to easily convert them into PAM values.
For 3-MF, our global BELGI-C s fit includes 19 parameters (see Table S2 of the supplementary materials): like for the 2-MF fit, we performed a fit on A, B, C, D ab , the five quartic and one sextic (H K ) centrifugal distortion constants, the parameter ρ, and the two first Fourier expansion coefficients V 3 and V 6 of the potential barrier. For 3-MF, we were able to also fit the F parameter and four higher order interaction terms between internal and global rotation. Our RMS deviations (181 kHz and 174 kHz for v t = 0 and 1 respectively) are slightly above the averaged measurement accuracy estimated to 150 kHz. Near 4500 rotational transitions involving energy levels up to J = 50 have been assigned in the 70-220 GHz scanned range. The 30 lines with low J values (18 A species lines and 12 E species lines) from Ogata and Kozima [11] with a weight of 100 kHz were also included to the fit. For those lines, the RMS deviation was 116 kHz. Sixteen lines from [11] showed large observed-calculated values of more than three times the measurement accuracy and we decided to discard them. As indicated in Table 2, for 3-MF the XIAM fit was of rather poor quality, with RMSs of 200.6 kHz and 1500 kHz for v t = 0 and v t = 1 respectively. However our BELGI-C s values for A, B, C transformed in the PAM only differ by 0.005, 0.04, 0.02% with those of XIAM in v t = 0. The ∠(i, a) angle differs by 0.2%. The V 3 values differs by 6% and 8% with XIAM v t = 0 and v t = 1, respectively. Again, as for the 2-MF, our PAM values fitted with the BELGI-C s code agree reasonably well with the B3LYP theoretical calculations.
For both the 2-and 3-MF, the threefold torsional potential barrier is rather high (V 3 ≈ 420 cm −1 and 358 cm −1 respectively), and the value of the reduced height s = 4V 3 /(9F) (with the internal rotation constant F ≈ 5.64 cm −1 for 2-MF and F ≈ 5.32 cm −1 for 3-MF) is 33 and 30, respectively. Torsional splittings reach values of 15 MHz in v t = 0 and 380 MHz in v t = 1 states. This classifies them as rather good candidates to apply "local" approaches at least in the v t = 0 state. Using the global approach, small internal rotation splittings may result into correlation between the parameters. As it was already observed in other small organic molecules (such as for methyl formate HCOOCH 3 [29]), the assignment was complicated by two main factors: (i) 2-MF and 3-MF have two small-amplitude vibrations [12]: the out-of-plane bendings of the methyl group at 237 and 267 cm −1 , respectively and the in-plane bending mode at 339 and 321 cm −1 , respectively. These modes are lying below the barrier and therefore are expected to cause some perturbations on the rotationtorsional energy levels. (ii) The value of the ρ parameter which represents the coupling term between internal rotation and global rotation in the kinetic energy operator is only 0.056. This low value does not allow us to get much information about the torsion-rotation coupling. In this complicated and rich room-temperature spectra of internal rotors, it is thus very useful to have various tools to first analyze the spectra (using the "local" codes which are very efficient and easy to use), and then fit globally the torsional manifold (with the "global" codes which are very powerful in predicting higher torsional states that can be populated). It is clear that so far no automatic procedure can really handle the complicated torsional-rotational manifold, and even more the vibrational-rotational-torsion manifold. A new code is currently built by Ilyushin and collaborators towards that goal [30].

Influence of the Isomerism on the Internal Rotation of the Methyl Group, Comparison with Other Heterocyclic Compounds
The hindered internal rotation of the methyl group in 2-MF and 3-MF is due to the one-dimensional torsional potential. The first two terms of its threefold Fourier series, see Equation (5), introduce a V 3 and a V 6 parameter. Quantum chemistry calculations based on Density Functional Theory (DFT) B3LYP method and the single reference Møller-Plesset perturbation theory at second order (MP2) method give the dependence of the energy of the electronic ground state on the torsional angle α. Figure 4 indicates that B3LYP and MP2 methods provide nearly identical barriers with slightly higher value with B3LYP. Furthermore, the barrier is lower in 3-MF than in 2-MF. This lowering of the barrier might be interpreted by the larger distance between the oxygen atom and the methyl group which slightly flattens the torsional potential. While the fits with just the V 3 parameter already qualitatively reproduce the data, the fits are improved by adding the V 6 term in the series expansion, in particular near the top of the barrier as it can be seen in the insets of Figure 4a,b and confirmed by the differences between the computed electronic energies and the fitted ones in Figure 4c,d.
The values of V 3 and V 6 obtained from quantum chemical calculations can be compared with the spectroscopic parameters of Table 1 for 2-MF and Table 2     fitting the experimental line positions. We observe that the V 3 parameter determined from quantum computations underestimates the XIAM or BELGI-Cs values by 10% for 2-MF (3% for 3-MF) while the V 6 parameter (in absolute value) underestimates the fitted values by a factor 5 or more for 2-MF (overestimates by a factor 4 for 3-MF). Such a disagreement is not surprising, because the V 3 and V 6 parameters obtained from the XIAM and BELGI-C s codes consider the averaging over the harmonic motion of the nuclei belonging to the furan frame and of the zero-point energy in the torsional potential. Because the previous measurements on 2-MF and 3-MF were focused on the ground state [11] or treat the v t = 0 and v t = 1 states independently in a local approach [12], those prior experimental studies dedicated to large amplitude internal rotation were only able to fit the V 3 barrier height without the higher-order V 6 correction. That was also the case for other heterocyclic analogs of methylfuran (MF) as methylthiophene (MT) and methylpyrrole (MP). In Table 3, we notice that the V 3 parameters were experimentally determined for the 2-and 3-isomers of each heterocyclic compounds but the higher-order V 6 correction could not be determined before our study of 2-MF and 3-MF mm-wave spectra since it requires the access to the v t = 1 excited torsional state and a global fit including both v t = 0 and v t = 1 allowed by the BELGI-C s code. Therefore no V 6 experimental values exist for the other heterocyclic compounds and the comparison should be performed with V 3 and V 6 values fitted on the theoretical data (hereafter B3LYP/aug-cc-pVTZ calculations). Table 3. V 3 and V 6 barriers for 2-and 3-methylfuran (MF), methylthiophene (MT) and methylpyrrole (MP) calculated with B3LYP/aug-cc-pVTZ as the function of the methyl groups torsional angle α. Uncertainties are the standard error given by Gnuplot [24] for the fit of the DFT potential based on Equation (5). The experimental values determined in this work for MF or in the literature for MT and MP are mentioned in square brackets. a Experimental value obtained from BELGI-C s . b Experimental value taken from [31]. c Experimental value taken from [32]. d Experimental value taken from [33]. e Experimental value taken from [34].

Molecule
The internal rotation potential and particularly the V 3 parameter is sensitive to both steric and electronic effects. In heterocyclic compounds, electronic effects are usually dominant [32,34]. When we compare the V 3 values of the three heterocyclic compounds, we can observe that: first, the variation of V 3 is stronger when the methyl group is substituted to the second position (2-MF, 2-MT and 2-MP) compared to its third position (3-MF, 3-MT and 3-MP); second, V 3 is larger for MF compared to MT and MP (( for the 3isomers) since a larger electronegativity of the oxygen which creates an electronic environment hindering the internal methyl rotation.
Concerning the V 6 term listed in Table 3, we observe larger values for the 3-isomers compared to the 2-ones as if the electronic effects tend to reduce the anharmonicity of the potential. Moreover we can observe negative V 6 values for the three 3-isomers where |V 6 /V 3 | > 3%. These negative V 6 are generally observed for internal rotation potential with significant anharmonicity, we can cite as example molecules such as acetaldehyde [35], methyl ketene [36] or meta-cresol [37].

Quantum Chemistry Calculations
The quantum chemical investigations on 2-MF and 3-MF were performed with the Gaussian 16 package [38]. The B3LYP [39] functional was used in the DFT calculations with the default ultrafine grid. We also performed molecular orbital based calculations with MP2 [40]. All our results were obtained with the correlation consistent basis set aug-cc-pVTZ [41]. The convergence criteria of the geometry optimisations of both 2-MF and 3-MF isomers were set to the tight option. The one dimensional energy curves describing the electronic energy as a function of the internal rotation of the methyl group are presented on Figure 4 and all results are listed in Table 3. Each energy value is the result of a geometry optimization on all the coordinates except the frozen angle describing the internal rotation. A step size of 5 • was chosen between two successive angles.

Internal Rotation Hamiltonian
There exist several approaches and codes to treat rotational spectra of molecules containing one or two internal rotor(s). Codes are available at the "Programs for ROtational SPEctroscopy" (PROSPE) [42], managed by Z. Kisiel. Among them, the three codes we are using in the present work, i.e., the program XIAM written by Hartwig and Dreizler [21] which can handle up to three internal rotors, the program called the "ERHAM" (ERHAM) developed by P. Groner for one or two-tops [22,27] and the BELGI series of codes developed by Hougen and Kleiner for asymmetric tops containing one or two internal rotors, with C s [43] (BELGI-C s ) and C 1 (BELGI-C 1 ) [29] symmetries. Finally, Ilyushin has developed the RAM36 code for one-top in the same RAM approach [44].
The Hamiltonian used in the BELGI-C s code is the RAM internal-rotation Hamiltonian based on the work of Kirtman [45], Lees and Baker [46], and Herbst et al. [47]. XIAM and ERHAM can be qualified as "combined axis methods" (CAM) because they set up the Hamiltonian (or parts of it) in the Rho-Axis-System (RAS) but then apply a transformation back to the PAS. ERHAM and XIAM allow to fit each torsional state separately (neglecting thus matrix elements off-diagonal in the torsional quantum number v t ). We thus classify these two methods as "local" ones because they fit each torsional state by itself, even though they treat the A and E states together. In the BELGI code, all the torsional states, up to a given truncation level are taken into account and fit together, so we call it a "global" method. This global approach is particularly successful for low barriers or high excited torsional states.
Since rather complete descriptions of these methods already exist in the literature [48], we will not repeat here such a general description. We only emphasize here the various characteristics that the three methods and programs had to face for the present study.
The RAM is a choice of an axis system which eliminates in the rotation-torsion Hamiltonian the Coriolis type terms (the x, y, z axis are related to b, c, a for prolate molecules. For an asymmetric top with a plan of symmetry as the case here, the Coriolis cross term −2Fρ y p α J y does not exist because of symmetry restrictions.) −2Fρ x p a J x and −2Fρ y p α J y , where J x and J y designate the components of the rotational angular momentum in the molecular axis system, p α is the internal rotational angular momentum, F is the internal rotation constant and ρ is the coupling constant. The RAM Hamiltonian may be written as [47]: where H T is the torsional Hamiltonian, H R the rotational Hamiltonian, H cd the usual centrifugal distortion Hamiltonian, and H int contain higher order torsional-rotational interaction terms: For this prolate asymmetric top, the reduction A and the representation I r were chosen. By grouping the terms and if x, y, z → b, c, a we get: The potential function which hinders the internal rotation of the methyl group is described by a threefold Fourier series depending on the torsional angle α: The relation between the rotational constants A, B, C in the principal axis system and the constants in the Rho Axis system is given by the angle θ RAM which turns the PAS into the RAS: From a fit of the spectrum, the following parameters are derived from the BELGI-C s code: V 3 , F, A RAM , B RAM , C RAM , D ab and ρ (plus higher-order parameters).
In the ERHAM code [22,27], the overall rotation-internal rotation Hamiltonian is set up in a product basis of symmetric rotor eigenfunctions with eigenfunctions of the internal rotation Hamiltonian set up in a Rho Axis system. It is assumed that the internal rotation Hamiltonian has already been solved. In the matrix elements of the effective Hamiltonian, the integrals over the internal rotation variables are expressed as Fourier series. Their coefficients are called "tunneling coefficients", and they are the fitting parameters of the effective rotational Hamiltonian in ERHAM.
In the XIAM code [21], the Hamiltonian is set up in the principal axes system of the entire molecule but then the internal rotation operator of each top is set up in its own Rho Axis system (RAS) and after diagonalization, the resulting eigenvalues are transformed (rotated) back into the principal axis system. The fitted parameters are the rotational constants A, B, C, the centrifugal distortion constants up to the 6th order, and for each methyl top the direction of the top axis in the PAM system ∠(i, x), ∠(i, y), ∠(i, z) the internal rotation constant F and the potential barrier height V 3 . Only some higher (4th order) coupling terms between internal rotation and overall rotation are sometimes implemented, and each torsional state is fit by itself.

Room Temperature Millimeter-Wave Spectroscopy
The Doppler limited room temperature mm-wave absorption spectra of 2-and 3-MF in the range of 70-110 GHz and 170-220 GHz were recorded at the LPCA laboratory using the solid state THz spectrometer described in detail in Ref. [20]. The measurements were performed in flux conditions with a low pressure of about P = 5 Pa allowing to limit the consumption of MF and to reduce the collisional broadening up to the Doppler broadening limit (∆ν Doppler ∼ 150 kHz (FWHM) at ν 100 GHz and T 300 K). Due to the high volatility of MFs (P vap 200 hPa for 2-MF and P vap 100 hPa for 3-MF at 300 K), the samples of 2-MF and 3-MF were cooled with a mixture of ice and liquid nitrogen and the flows were regulated with a micro-valve. The spectra were recorded using a combination of a frequency multiplier chain (Virginia Diodes, Inc., Charlottesville, VA, USA) exciting the molecule and a ZBD, which was an unbiased Schottky diode mounted in a wave-guide operating in detection mode.
After up-converting the synthesized microwave frequency by the multiplication chain the guided radiation is launched into free space using a horn antenna, collimated by an off-axis parabolic mirror, passing a polarization grid and propagated through a 125 cm long and 5.6 cm diameter stainless steel absorption cell closed by two Teflon windows. The interaction path-length is doubled [49] with a roof-top reflector positioned at the output of the gas cell. The reflected beam is returned to the output parabolic mirror which ensures the focusing onto the detector. The second harmonic of the frequency modulated signal were recorded in 100 kHz steps with a time constant of 100 ms using a modulation depth of 50 kHz and a modulation frequency of 43 kHz delivered by the frequency synthesizer.
A post treatment with a FFT band-pass filter cuts off lower frequencies and the noisy higher frequency parts of the spectra reducing the noise level and removing the low frequency baseline variations caused by stationary waves between the source, the polarization grid, the cell, and the detector [28].

Conclusions
2-MF and 3-MF isomers are oxygenated VOCs emitted from multiple sources in the atmosphere especially natural or industrial biomass combustion processes. Their atmospheric reactivity leads to the formation of SOAs which directly or indirectly affect the climate. Mm-wave spectroscopy is an interesting approach to monitor them but it requires a prior determination of reliable rotational-torsional line lists. In this study, we propose a global approach using the BELGI-C s code providing a unique set of fitted parameters able to reproduce, at the experimental accuracy, the 2-MF and 3-MF observed mm-wave lines belonging to the ground v t = 0 and the first torsional v t = 1 states up to high J and K values. Unlike 2-MF for which our assignment was based on previous measurements performed by Finneran et al. in Ref. [12], the mm-wave spectrum of 3-MF has been measured in this study in the 70-220 GHz frequency range with the mm-wave spectrometer developed at the LPCA [20]. For both isomers we compared our global approach using BELGI-C s with local approaches based on the XIAM and ERHAM codes. Two advantages of the BELGI-C s global approach are highlighted: (i) the ability to fit at the experimental accuracy the A and E components involving both v t = 0 and v t = 1 states (XIAM failed to fit correctly rotational lines of 2-MF and 3-MF in v t = 1); (ii) the possibility to determine with a high degree of accuracy the V 3 and V 6 parameters of the internal rotation potential (in the ERHAM method, the potential barrier does not appear explicitly and V 6 could not be determined). Based on DFT calculations, these (V 3 , V 6 ) parameters are discussed considering other sulfur and nitrogen heterocyclic analog in order to better understand the influence of the electronic environment on the internal rotation barrier and its higher-order correction term. To confirm experimentally the identified trends, we propose to investigate in the mm-wave domain the rotational-torsional spectroscopy of MP and MT already studied in the microwave domain [31][32][33][34]. Finally, in order to assign the remaining mm-wave lines of 2-MF and 3-MF, a rovibrational analysis of the lowest energy vibrational modes will be undertaken by means of synchrotron-based FT-Far-IR high-resolution spectroscopy as it was done recently by our group on the catechol molecule [15].