Toward Accurate yet Effective Computations of Rotational Spectroscopy Parameters for Biomolecule Building Blocks

The interplay of high-resolution rotational spectroscopy and quantum-chemical computations plays an invaluable role in the investigation of biomolecule building blocks in the gas phase. However, quantum-chemical methods suffer from unfavorable scaling with the dimension of the system under consideration. While a complete characterization of flexible systems requires an elaborate multi-step strategy, in this work, we demonstrate that the accuracy obtained by quantum-chemical composite approaches in the prediction of rotational spectroscopy parameters can be approached by a model based on density functional theory. Glycine and serine are employed to demonstrate that, despite its limited cost, such a model is able to predict rotational constants with an accuracy of 0.3% or better, thus paving the way toward the accurate characterization of larger flexible building blocks of biomolecules.


Introduction
The experimental study of biomolecule building blocks in the gas phase has recently attracted increasing attention owing to the development of spectrometers coupling supersonic-jet expansion [1] with laser ablation [2]. This has allowed the application of high-resolution rotational spectroscopy in the microwave region to thermolabile molecules with high melting points. Among these, we find the main families of life bricks such as amino acids, peptides, DNA and RNA bases. The focus of the present contribution is the computational support to rotational spectroscopy investigations of building blocks of biomolecules.
For flexible molecules such as those of interest for this study, a great challenge is related to the large number of conformers and the fast relaxation of some of them to more stable counterparts due to the presence of low interconversion energy barriers. An inaccurate account of the relaxation processes can bias any direct thermochemical interpretation of the results provided by rotational spectroscopy experiments [3,4]. Quantum-chemical (QC) computations can help to tackle this challenge, especially because gas phase is their most natural playground [5,6]. Unfortunately, already for medium-sized systems, the usual dichotomy between accuracy and feasibility, which is the quest for accurate yet feasible predictions, comes into place [7]. In fact, state-of-the-art QC approaches are able to rival the experimental counterparts for small semi-rigid systems in the gas phase [7][8][9], but they are characterized by a very unfavorable scaling with the dimension of the system to be investigated. This prevents their brute-force application already to biomolecule building blocks containing more than a dozen atoms and characterized by several low-energy minima. Furthermore, the powerful local optimization techniques developed for semi-rigid systems are not effective for flexible systems, which requires the exploration of rugged potential energy surfaces (PESs) [10,11].
For the reasons mentioned above, the accurate characterization needed by rotational spectroscopy requires an integrated computational approach that employs QC models of increasing accuracy in the different steps of an exploration/exploitation strategy guided by machine learning (ML) tools. The main steps of this strategy [10,[12][13][14] can be summarized as follows: 1.
Unsupervised perception of the molecular system to identify hard and soft degrees of freedom [15]; 2.
Knowledge-based selection and constrained geometry optimizations of a limited number of conformers employing a fast semi-empirical method [11,16]; 3.
Exploration of the PES governed by soft degrees of freedom using the same semiempirical method of the previous step, guided by a purposely tailored evolutionary algorithm with the aim of finding other low-lying minima [10]; 4.
Refinement of the most stable structures by hybrid and then double-hybrid density functionals [14]; 5.
Analysis of relaxation paths between pairs of adjacent energy minima [13]; 6.
Computation of spectroscopic parameters for the energy minima with non negligible populations [13,26].
The focus of the present contribution is on points 6 to 8. First, we will validate the accuracy of approximate approaches based on double-hybrid functionals and effective composite methods by exploiting results available for semi-rigid molecules and non-covalent complexes. Next, amino acids will be employed to illustrate the potentialities of the validated approaches. These compounds represent a particularly appealing playground because their rich conformational landscape is tuned by the competition among different types of intra-molecular non-covalent interactions involving, together with the amino and carboxylic acid moieties of the backbone, also side-chain groups. At the same time, results from MW experiments are available for several conformers of most natural α-amino acids [27][28][29][30][31][32][33][34][35][36][37][38], and provide accurate data for benchmarking theory. Among the α-amino acids well characterized experimentally, we have selected glycine and serine.
Since the current standard for the computation of rotational spectroscopy parameters of biomolecule building blocks (see, for example, [4,27,[39][40][41]) is based on QC methodologies of limited accuracy and does not account for vibrational effects, together with the intrinsic interest of the studied molecules, the results of the present study will provide a reference for a more accurate and reliable approach to be used in the prediction and analysis of MW experiments.

The Methodologic Approach
Several studies have shown that the double-hybrid rev-DSD-PBEP86 functional [42] in conjunction with the jun-cc-pVTZ [43] basis set (this functional-basis combination is shortly denoted as rDSD in the following) represents a very robust and convenient choice for obtaining equilibrium geometries suitable for subsequent electronic energy evaluations [19,[44][45][46]. For the latter, composite schemes rooted in the coupled-cluster (CC) theory [47] that account for the extrapolation to the complete basis set (CBS) limit and core-valence (CV) correlation are employed whenever aiming at high accuracy [13,25,46,[48][49][50]. However, composited schemes entirely based on CC theory are computationally expensive. The key idea of the so-called 'cheap'-scheme (ChS) family of composite approaches [17][18][19]51] is that both CBS and CV contributions can be estimated accurately by low-order (hence low-cost) many-body perturbative methods, with Møller-Plesset theory to second order (MP2) [52] offering the best compromise between accuracy and computational cost. The general expression of the electronic energy evaluated using the ChS model is: where the first term on the right-hand side is the energy at the fc-CCSD(T) [53] level in conjunction with a triple-zeta (TZ) basis set, where CCSD(T) denotes the CC method incorporating single, double, and a perturbative estimate of triple excitations, and fc stands for the frozen-core approximation. The second term is the contribution due to the extrapolation to the CBS limit evaluated at the MP2 level within the fc approximation and using the n −3 formula [54]: where n = 4 (quadruple-zeta basis set, QZ) and n-1 stands for TZ. The last term, ∆E(MP2/CV), is the CV contribution computed as energy difference between MP2 calculations correlating all electrons and within the fc approximation, both in the cc-p(w)CVTZ basis set [55,56]. The first implementation of the ChS model employed the cc-pVnZ basis sets [57]. However, several benchmarks [18,19] have led to the conclusion that the best cost/performance ratio is obtained with the jun-cc-pVnZ basis sets [43], with jun-cc-pV(n+d)Z being used for atoms of the third period [58]. This model, hereafter referred to as junChS, has become our standard ChS methodology. A more recent alternative [19] replaces the conventional CCSD(T) and MP2 methods with their explicitly-correlated counterparts, namely CCSD(T)-F12 and MP2-F12 [59,60], still in conjunction with the jun-cc-pVnZ and cc-p(w)CVTZ basis sets. Hereafter, this model is referred to as junChSF12. In all CCSD(T)-F12 computations, the F12b approximation [61] is employed, while the default approximation (3C/FiX) is used for MP2-F12 [60].
In addition to energy evaluations, different ChS variants have been employed and tested for geometry optimizations [19,51,62,63]. For conventional QC methods, ChS and junChS are defined analogous to Equations (1) and (2), whereas the extrapolation to the CBS limit (but not the CV correction) can be avoided (and sometimes introduces nonnegligible oscillations) for approaches exploiting explicitly-correlated methods. The models derived and employed in the present work are denoted as DZCCF12+CV, TZCCF12+CV, junCCF12+CV and augCCF12+CV, which correspond to CCSD(T)-F12 calculations in conjunction with cc-pVDZ-F12, cc-pVTZ-F12, jun-cc-pVTZ and aug-cc-pVTZ, respectively. The "+CV" denotes the incorporation of the CV contribution, which is evaluated at the MP2-F12/cc-p(w)CVTZ level.
When dealing with flexible molecules, a key point is the determination of the relative stability of the low-energy minima to determine those that are sufficiently populated for a spectroscopic characterization. To accomplish this task, one has to move from electronic energy differences to the corresponding relative enthalpies or free energies, evaluated at a temperature that depends on the experimental conditions. To this end, ZPEs and vibrational partition functions are needed. Within the harmonic approximation, this task requires the evaluation of harmonic vibrational frequencies, which can be obtained with a suitable accuracy at the rDSD level [64] exploiting analytical second derivatives [65]. The models of the ChS family can also be applied to the computation of harmonic frequencies.
However, due to effective error compensation, improved results are obtained by neglecting the CV correction for explicitly-correlated approaches [66]. Focusing on models that exploit the F12 methodology, the DZCCF12, TZCCF12, junCCF12 and augCCF12 approaches offer a good compromise between accuracy and computational cost. The ZPE term can then be improved by incorporating anharmonic corrections, which are usually obtained at the B3LYP/jun-cc-pVDZ level (hereafter B3) within second-order vibrational perturbation theory (VPT2, [20,67]). Whenever required, the contribution of low-frequency motions to en-tropies can be computed using the quasi-harmonic (QH) approximation [24,68]. Within the QH approximation, below a given cut-off value, entropic terms are obtained from the free-rotor approximation, and a damping function is used to interpolate between free-rotor and harmonic oscillator expressions close to the cut-off frequency.
Moving to rotational spectroscopy, the leading terms are the rotational constants, which incorporate two contributions. The first, and by far the largest, one is the equilibrium rotational constant, B e , which is straightforwardly derived from the equilibrium geometry of the studied molecule. The second contribution is related to vibrational effects which are present because of the coupling between vibrations and rotation. This contribution is usually approximated by means of VPT2, which, for the vibrational gound state, leads to: where α i r denotes the vibrational-rotation interaction constants, with i being the inertial axis (a, b, or c) and the sum running over all the r vibrational modes. The computation of this contribution requires both second-and semi-diagonal third-energy derivatives of the energy with respect to normal modes. Therefore, computations of vibrational corrections to rotational constants are quite expensive. As a consequence, its inclusion is warranted only if the error on the computed equilibrium rotational constants is smaller than the expected vibrational contribution. A recent careful analysis of these aspects provides some general hints [69]. The magnitude of the vibrational contribution ∆B vib is typically 0.1% to 0.7% that of the corresponding equilibrium rotational constant, with 0.5% being a very robust guess for semi-rigid molecules. As a consequence, an uncertainty of 10% on ∆B vib would lead to an error of 0.05% on the rotational constants, which is thus more than acceptable. The overall conclusion is the suitability of global hybrid density functionals in conjunction with double-zeta basis sets (e.g., the B3 model previously defined) for the computation of the required anharmonic force fields, thus largely reducing the computational cost.
Focusing on the accuracy needed by rotational spectroscopy, errors of 1% on rotational constants (100 MHz for a constant with magnitude 10 GHz) are not at all helpful for the prediction and/or analysis of rotational spectra of flexible molecules, which are always characterized by several low-energy minima with non negligible populations. On the other hand, an accuracy of 0.01%, which would allow direct comparison with experiment, is presently reachable only for molecules containing two or three atoms. Based on these considerations, the optimal level of accuracy associated with predicted rotational constants should be close to 0.1% (10 MHz for a constant of 10 GHz). In terms of structural parameters, such an accuracy corresponds to errors smaller than 0.001 Å for typical bond lengths and 0.001 radians (0.05 degrees) for typical valence angles [69]. This target accuracy can be surely obtained by expensive composite schemes incorporating high excitation orders in the correlation treatment [70]. Nevertheless, the models of the ChS family are able to draw closer to such a high accuracy limit [51,62,69]. However, for quite large flexible molecules, even this computational level becomes too expensive and one has to rely on density functional theory (DFT).
To improve the optimized geometries obtained from double-hybrid functionals, one can resort to the the so-called linear regression approach (LRA) [71][72][73][74]. This corrects the computed bond lengths (r comp ) for systematic errors by means of scaling factors (a) and offset values (b) that have been derived from the comparison of DFT and accurate semi-experimental (SE) equilibrium geometries for a large database: with the a and b parameters depending not only on the functional considered, but also on the nature of the atoms involved. These are available in [73]. As mentioned above, the equilibrium rotational constants contribute more than 99% to the parameters derived from experiments and only depend on the molecular structure and isotopic composition [75]. Indeed, rotational constants are inversely proportional to the inertia tensor, which contains information on the mass distribution in the molecule. Therefore, a comparison between computed and experimental rotational constants can provide hints on the quality of the computed structure. In this respect, additional parameters of particular relevance are the nuclear quadrupole coupling constants (χ ii , i referring to the inertia axis a, b or c), which strongly depend on the intra-molecular interactions [76]. Nuclear quadrupole coupling is the interaction between the quadrupole moment of a quadrupolar nucleus and the electric gradient at the nucleus itself, with quadrupolar nuclei being those having a nuclear spin I ≥ 1 [75]. Since 14 N is a quadrupolar nucleus and is present in almost all biomolecule building blocks, nuclear quadrupole coupling constants can be used, together with rotational constants, to assess the quality of computed structures. In turn, their calculation is important for accurate predictions of rotational spectra because nuclear quadrupole coupling determines a splitting of the rotational energy levels, thus causing a splitting of the rotational transitions, which is the so-called hypefine structure. It is also noted that vibrational effects on nuclear quadrupole coupling constants are usually smaller than the uncertainty affecting the computed equilibrium values, and thus have not been considered in this work.
As far as technical details are concerned, all DFT, MP2 and single-point energy CCSD(T) computations have been performed using the Gaussian program package [77], while the CFOUR program [78] has been employed for geometry optimizations at the CCSD(T) level. For explicitly-correlated calculations, the MOLPRO program [79] has been used.

The Validation Step
Non-covalent complexes represent a suitable benchmark for the validation of computational methodologies aiming at an accurate description of intra-molecular non-covalent interactions. In a series of recent papers [18,19,80], different ChS variants have been applied to the computation of the interaction energies of prototypical systems. These studies (partially summarized in Table 1) show that the junChS and, especially, junChSF12 approaches allow the determination of interaction energies with an average error smaller than 10 cm −1 (0.12 kJ mol −1 ) without the need for incorporating any empirical parameter. In Table 1, the non-covalent complexes of the A14 database [19], which is a based on a selection of systems from the A24 dataset [81], have been considered. As reference values, the highly accurate data (which also incorporate up to quadruple excitation in the CC expansion and relativistic effects) from [82] have been selected.
From an inspection of Table 1, we note that going from junCCF12+CV to junChSF12, which only differ for the extrapolation to the CBS limit, leads to a significant lowering of both the maximum and mean unsigned error (MAX and MUE, respectively), with this confirming the importance of the extrapolation to the CBS limit also for explicitly-correlated methods and that this contribution can be accurately incorporated by means of low-order perturbative methods (here MP2-F12). Indeed, MP2-F12 represents a reliable route for obtaining accurate results without any significant computational increase with respect to the underlying CCSD(T)-F12 step.
Despite the extrapolation to the CBS limit, all ChS models (both conventional and F12 variants) are affected by small, but not entirely negligible basis set superposition errors (BSSE), which can be taken into account by means of counterpoise (CP) corrections (for details, the reader is referred to [19]). However, these are slightly more pronounced when considering the model based on conventional methods. This is particularly significant for intra-molecular interactions where the BSSE becomes ill defined and CP corrections difficult to implement. An important outcome of [18,19] is that, when using the rDSD level for the evaluation of the reference structures, non-CP corrected geometries can be safely employed. This is particularly important in view of applying junChS and junChSF12 to the evaluation of the energetics for flexible molecules. Concerning interaction energy, the errors at the rDSD level without any CP correction nearly double with respect to junChS(F12) [19], but they remain significantly smaller than those issued by other density functionals. More importantly, the rDSD level is suitable for defining a correct stability order for different kinds of non-covalent interactions. a Energies evaluated at the reference geometries from [81]. b Errors are with respect to the best values from [82]. c MAX, MUE and RMSD stand for maximum error, mean unsigned error, and root-mean-square deviation, respectively.
As mentioned above, the largest contribution to rotational constants comes from their equilibrium value, which is related to the corresponding equilibrium geometry. In this connection, ChS, junChS and junChSF12 have proven to provide very accurate results [19,63], with the former two models being rather well tested [83]. When employing explicitly-correlated methods, the computational cost of the ChS methodology can be further reduced by omitting the CBS extrapolation step. To inspect this, a validation of the junChS and CCF12+CV models (for CCF12, we considered the DZCCF12, TZCCF12 and junCCF12 variants) has been performed by comparing their results with accurate SE equilibrium geometries from the SE100 database [73]. The subset considered includes 24 molecules containing hydrogen, second-and third-row atoms: CH 4 , CO 2 , HCN, HNC, H 2 O, NH 3 , C 2 H 2 , C 2 H 4 , H 2 CO, t-HCOOH, CH 2 NH, BH 3 NH 3 , BH 2 OH, C 2 H 4 O (oxirane), C 2 H 4 NH (aziridine), C 3 H 6 (cyclopropane), H 2 O 2 , SO 2 , H 2 S, PH 3 , H 2 CS, CH 2 PH, C 2 H 4 S (tiirane) and H 2 S 2 . The results, collected in Table 2 in terms of maximum and mean unsigned errors (MAX, MUE), show that the CV contribution cannot be neglected and that the jun-cc-pVTZ basis set is competitive with its counterparts (cc-pVDZF12 and cc-pVTZF12) developed specifically for explicitly correlated models.
Even if not reported in Table 2, from [73], it is observed that-for the entire SE100 setthe rDSD model delivers remarkably accurate geometrical parameters with typical absolute errors of about 0.003 Å for bond lengths and 0.1-0.2 degrees for valence angles. Furthermore, the errors appear to be very systematic. Therefore, the rDSD level is particularly suitable for the application of the LRA mentioned above, which greatly reduces the expected errors [73]. In particular, in view of the focus of this work, the accuracy for the corresponding rotational constants becomes close to the 0.1% uncertainty-target.  [90] and PH 3 [91]). As already pointed out, CV contributions have not been considered.
As before, results are reported in terms of MAX, MUE and RMSD. From Table 3, it is quite apparent that the different basis sets tested in conjunction with the CCSD(T)-F12 ansatz provide very similar results. Therefore, the jun-cc-pVTZ basis set, which has been found to be the optimal choice for energy and geometry evaluations, can be confidently employed also for harmonic frequencies. Furthermore, the much cheaper rDSD level performs remarkably well and can be profitably used for the computation of harmonic frequencies of larger molecules and thus for the accurate evaluation of their ZPEs.

Amino Acids
The structure of isolated amino acids is ruled by both backbone (φ = HNC α C', ψ = NC α C'O(H) and ω = C α C'OH) and side-chain (χ) torsional angles (see Figure 1 for the specific case of serine). The non-planarity of the NH 2 moiety suggests that, instead of the customary φ dihedral angle, φ = LP-N-C α -C' = φ + 120 • can be better employed. In the definition above, LP is the nitrogen lone-pair perpendicular to the plane defined by the two aminic hydrogens and the C α atom. The most stable backbone structures result from the formation of hydrogen bonds, which can be classified as follows, also according to the values of the dihedral angles defined above. Type I is a bifurcated NH 2 · · · O=C hydrogen bond (φ ≈ 180 • , ψ ≈ 180 • , ω ≈ 180 • ), type II is a N· · · HO hydrogen bond (φ ≈ 0 • , ψ ≈ 0 • , ω ≈ 0 • ), and type III is a bifurcated NH 2 · · · OH hydrogen bond (φ ≈ 180 • , ψ ≈ 0 • , ω ≈ 180 • ). The type of hydrogen bond in the backbone also leads to the principal denomination of the structure, as shown in Figure 2. It has to be noted that type III structures have never been observed in MW experiments of natural α-aminoacids without polar side-chains due to their easy relaxation to the corresponding type I counterparts. Additional hydrogen bonds can be established with polar side chains, and they further stabilize the structure. This is the case, for example, of the type III' conformers (single NH· · · OH hydrogen bond, φ ≈ 180 • , ψ ≈ 90 • , ω ≈ 180 • ). Typical conformers of types I, II and III' are shown in Figure 2 for the case of serine (which will be discussed in detail in the following), where the g − , g and t labels refer to the gauche and trans conformations of each χ dihedral angle (as evident from Figure 1, for serine, there are two χ angles).
While the strategy described in the Introduction is able to identify all possible conformers of a given amino acid, their experimental characterization is limited by the sensitivity of the technique employed. For the specific case of rotational spectroscopy, a conservative limit for the relative stability of detectable structures is around 900 cm −1 , which corresponds to a relative population of about 1% at room temperature (where kT/hc = 207 cm −1 ). In the same vein, an upper limit of 400 cm −1 for interconversion barriers is usually employed for discriminating fast relaxation processes in the case of flexible compounds showing several stable minima [4,39].
For the simplest amino acid, i.e., glycine, eight conformers have been characterized computationally, with only two of them (I and II) being sufficiently stable (and non-involved in fast relaxation processes) to be detected in MW experiments [92]. The limited size of glycine has allowed the evaluation of relative stabilities by means of state-of-the-art composite schemes including, together with CBS and CV contributions evaluated at the CCSD(T) level, full account of triple excitations, perturbative estimate of quadruple excitations, diagonal corrections to the Born-Oppenheimer approximation and relativistic effects [93]. While the reader is referred to [49,[92][93][94] for details, here we point out the remarkable performance of the ChS and junChSF12 approaches in the evaluation of relative energies (Table 4), with MUEs smaller than 10 cm −1 with respect to the most accurate results [93].
Noticeably, the rDSD model performs well with a MUE of 15 cm −1 . It is worth noting that such an accuracy is by far better than that obtainable by the methods (B3LYP and MP2) usually employed in the interpretation of MW spectra of biomolecule building blocks (see Table 4). Similar remarks apply to spectroscopic parameters. In fact, the MUEs with respect to experimental rotational constants are 23.6 MHz for rDSD, 13.0 MHz for rDSD-LRA, and 17.1 MHz for ChS. These average deviations from experiment point out the great accuracy that can be reached by applying LRA corrections to rDSD structures. For glycine, the availability of an accurate SE equilibrium structure for conformer I permits to evaluate directly the accuracy of geometrical parameters for the rDSD and rDSD-LRA levels. For rDSD, the MUEs are 0.0019 Å for bond lengths and 0.0026 radians (0.15 degrees) for valence angles, with the former value decreasing to 0.0004 Å at the rDSD-LRA level. Moving to nuclear quadrupole coupling constants, we note the good agreement (nearly quantitative) with experiment obtained by the ChS and rDSD levels, with the latter requiring a greatly reduced computational effort with respect to the former. This outcome further confirms the good performance of the rDSD model. Table 4. Equilibrium rotational constants and 14 N-nuclear quadrupole coupling constants (MHz) of the two most stable conformers of glycine and serine. The relative energy differences (∆E, in cm −1 ) are also reported.  [49] for glycine and [95] for serine. Experimental nuclear quadrupole coupling constants are from [28] for glycine and [33] for serine. b Glycine: ChS and MP2 results from [49]; rDSD and rDSD-LRA results from this work; B3LYP results from [71] and this work. Serine: ChS results from [95]; all the other data from this work. Systematic investigations have revealed that, analogous with glycine, the natural amino acids containing simple non-polar side chains (alanine [29], valine [30], isoleucine [31] and leucine [32]) have two dominant conformers of type I and II. On the contrary, the conformational landscape of natural amino acids with polar side chains is much richer due to the cooperativity or competition between intra-backbone and backbone with side-chain hydrogen bonds. The simplest amino acid containing a polar side-chain is serine, which has two degrees of freedom in its CH 2 OH side-chain: χ 1 = N-C α -C β -O and χ 2 = C α −C β -O-H (see Figure 1). The exploration of the conformational PES and refinement of low-lying energy minima end up with 12 structures within 900 cm −1 above the absolute energy minimum. However, five of them relax to more stable counterparts through low energy barriers. Thus, seven conformers are left as possibly detectable in MW experiments: three of type II, two of type III', and one each for types I and I'. Contrary to the case of glycine, the most stable conformer of serine in terms of electronic energies is of type II. However, incorporation of ZPE and thermal effects reverses the situation. As a result, the most stable conformer is Ig − g in terms of standard free energies.
The rotational constants of the two most stable conformers have been recently computed by geometry optimizations at the ChS level, with relative MAX and MUE of 0.6% and 0.3% with respect to the experiment, respectively [95]. It is noteworthy that very similar (but slightly improved) relative errors (0.4% and 0.3%, respectively) are obtained at the rDSD-LRA level (see Table 4). The strongly reduced computational requirements of this level have allowed us to obtain accurate results also for the other five low-energy conformers experimentally characterized, the results being provided in Table 5. It is noted that the comparison with the experiment confirm (i) an average relative deviation smaller than 0.4% for rotational constants and (ii) the good accuracy that rDSD is able to reach in the prediction of quadrupole coupling constants. This table also points out that the first four conformers (one of type I, one of type I' and two of type II) are significantly more stable than the remaining ones, with the Ig − g conformer becoming the most stable one once moving from relative electronic energies to standard free energies. At the rDSD/B3 level (rDSD referring to the level of theory used for electronic energies, B3 to the level employed for the computation of ZPE and thermal corrections), the IIgg conformer lies 44.1 cm −1 above (38.7 cm −1 employing ChS electronic energies). Table 5. Ground-state rotational constants (A 0 , B 0 , and C 0 in MHz) and 14 N-nuclear quadrupole coupling constants (χ ii , i = a, b, c, in MHz) of the seven conformers of serine experimentally investigated. Computed relative standard free energies (in cm −1 ) are also reported.  (21) a All the computed data are at the rDSD level (rDSD-LRA for equilibrium rotational constants) except for vibrational corrections to equilibrium rotational constants, ZPE and thermal contributions to ∆G, which are all at the B3 level. b Experimental data are from [33]. Standard errors are shown in parentheses in units of the last digits.

Conclusions
In this paper, we have focused on the final steps of a general strategy that aims at the accurate structural and spectroscopic characterization of flexible molecules. The main outcome of this work is that accurate geometries, spectroscopic parameters and relative energies can be obtained by the double-hybrid rev-DSD-PBEP86 functional. In particular, when combined with the jun-cc-pVTZ basis set, this DFT model is able to provide results approaching the accuracy of those issued by the well-tested ChS approaches that exploit the gold-standard CCSD(T)(-F12) method. The performance of the ChS composite schemes for geometrical parameters (and thus equilibrium rotational constants) is matched (if not improved) by rDSD when this is coupled with the linear regression approach (LRA), which does not involve any additional computational cost.
In detail, we have demonstrated that rotational and quadrupole coupling constantsparameters that strongly depend on mass distribution and intra-molecular interactionscan be obtained at the rDSD level with an accuracy which is suitable for guiding and supporting experiments in the field of rotational spectroscopy. The results obtained for glycine and serine are in very good agreement with the available spectroscopic data and pave the way toward the accurate investigation of large biomolecule building blocks.