Modelling the non-local thermodynamic equilibrium spectra of silylene (SiH2)†
This paper sets out a robust methodology for modelling spectra of polyatomic molecules produced in reactive or dissociative environments, with vibrational populations outside local thermal equilibrium (LTE). The methodology is based on accurate, extensive ro-vibrational line lists containing transitions with high vibrational excitations and relies on the detailed ro-vibrational assignments. The developed methodology is applied to model non-LTE IR and visible spectra of silylene (SiH2) produced in a decomposition of disilane (Si2H6), a reaction of technological importance. Two approaches for non-LTE vibrational populations of the product SiH2 are introduced: a simplistic 1D approach based on the Harmonic approximation and a full 3D model incorporating accurate vibrational wavefunctions of SiH2 computed variationally with the TROVE (Theoretical ROVibrational Energy) program. We show how their non-LTE spectral signatures can be used to trace different reaction channels of molecular dissociations.
Introduction<!>Disilane isomers<!><!>Disilane isomers<!>The silylene fragments<!><!>Modelling the non-LTE populations of SiH2<!>Computing vibrational populations of SiH2<!>The 1D approach<!>1D populations and spectra<!>The 3D approach for vibrational populations using an accurate variational method<!>PESs of SiH2 fragments<!>Vibrational calculations<!>Vibrational populations of SiH2<!>Non-LTE intensity simulations<!>Using the non-LTE populations from the 3D approach<!>Ã–X̃ spectrum<!>Conclusion<!>Conflicts of interest
<p>Normally, molecules are assumed to be in local thermal equilibrium pertaining to a given temperature with the internal degrees of freedom (electronic–rotation–vibration) characterized by the Boltzmann distribution. However many different physical chemical, experimental and technological processes produce molecules that do not satisfy the Boltzmann law and as consequence have unusual, non-local thermal equilibrium (non-LTE) spectroscopic signatures. Molecules produced in reactions do not necessarily obey the Boltzmann thermal equilibrium, at least if the reaction time is shorter than the collision time. Instead, their internal degrees of freedom are populated based on the reaction paths rather than on the temperature of the surrounding environment. These out-of-LTE (i.e. non-LTE) populations encode information about the structural reaction dynamics and can manifest in the molecular spectra. The field of non-LTE spectroscopy has great potential to study these processes as the properties of the molecules producing the non-LTE spectroscopic signatures can shed the light on the dynamics of chemical reactions.1–6 The so-called transition state (TS) spectroscopy is a technique already widely used that employs the high-resolution non-LTE spectra of products to observe reaction processes that are hidden for the conventional spectroscopic methods.7,8 The novel high resolution non-LTE spectroscopic techniques allow decoupling of the vibrational and rotational degrees of freedom of molecules and thus control their vibrational and rotational populations, e.g. with rich vibrational and simplified rotational structures.6,9–13</p><p>The modern day study of non-LTE spectroscopy can be traced as far back as the 1930s to the original papers of Milne,14 and the many key papers from the decades following.15–24 The non-LTE spectroscopic effects play important role in high-resolution applications and there exist a number of accurate non-LTE spectroscopic and radiative transfer codes, see van der Tak et al.,25 Funke et al.,26 Pannier and Laux27 and references therein. As such, non-LTE spectra are often vital for the modelling of astrophysical problems, including planetary atmospheric properties,28 stellar atmospheres of solar system and exoplanets29–36 and the ISM.25,37,38</p><p>The non-local thermodynamic effects of the spectra of molecules has been of interest to chemists and astronomers alike for many years. A notable example of this is the 2011 work by Ferus et al.2 who studied the isomers of HCN within acetonitrile, formamide, and BrCN discharge. Using the features for both the HCN and HNC molecules from these spectra, Ferus et al.2 were able to calculate the ratios of molecules within the reactions and also the reaction path taken by the HNC molecule during the isomerization. In this work we explore this idea to study non-LTE spectral signatures of silylene (SiH2) produced from disilane (Si2H6).</p><p>Reaction properties of silylene, silane and disilane such as the rate constants for the formation, destruction and chemical pathways are important for plasma physics aspects such as silicon deposition.39–41 The ease of hydrogen transfer and high barriers in the saturated silicon system, leading to the ready formation of three-centre interactions and consequently the isomerisation reactions of Si2H6, are just as important to study as the elementary reactions.42 Silane containing reactions are also of importance for astrophysics, with the presence of SiH4 in IRC + 10216 discussed by Goldhaber and Betz,43 Kaiser and Osamura.44</p><p>The complexity of the silane containing systems has been discussed elsewhere,42,45 with quantitative calculations proving particularly difficult. There were a number of ab initio studies of the structural properties of Si2H6e.g.,46–53 as well as of the formation and reactions involving this molecule.54–58 Agrawal et al.53 produced a global ab initio potential energy surface of disilane and used it to investigate the dissociation dynamics with classical trajectories. Márquez et al.50 reported a force field for Si2H6. The main isomer of Si2H6 has a staggered, ethane-like, structure (see Fig. 1) with a low barrier (1.2 kcal mol−1) to the eclipsed conformation.59 Si2H6 has been shown to have a local minimum as an inverted stable structure with one of the Si–H3 'umbrellas' pointing to the center as well as a transition state with a similar inverted configuration.46,60 These structures are nominally asymmetric (Cs symmetry) but essentially acquiring the C3v symmetry.</p><p>Thermal decomposition of Si2H6 has been extensively studied, both theoretically (mostly using RRKM, Rice–Ramsperger–Kassel–Marcus) and experimentally, with the reaction Si2H6 → SiH2 + SiH4 as the most common53,56,57,61 and important decomposition process of the excited disilane,53,57,62–66 and where Arrhenius parameters and rate constants have been reported (e.g. Bowrey and Purnell,61 Martin et al.,63 Mick et al.,65 Roenigk et al.67). However, it is also possible for the disilane molecule to dissociate homolytically, as was originally thought to be the main pathway owing to disilanes similarities with ethane, and form ,53,56,68 or to undergo dehydrogenation to H2Si–SiH2, H3Si–SiH, or H3Si–H2.53,54,57 Disilane can even undergo double dehydrogenation to form Si2H2, however to our knowledge this has only been reported as the main product when undergoing photolysis at 193 nm.69 Yoshida et al.57 notes that the transition state for Si2H6 → SiH2 + SiH4 is 8.48 kcal mol−1 lower than the transition state for Si2H6 → H3SiSiH + H2, at 43.38 kcal mol−1 compared to 51.86 kcal mol−1.</p><p>The spectroscopy of SiH2 has been used to monitor the SiH2 + SiH4 → Si2H6 and Si2H6 → SiH2 + SiH4 reactions (see Fig. 1) and measure the corresponding rate constants and Arrhenius parameters by spectroscopically tracking electronic (Ã–X̃) transitions of SiH2.39,46,65,70 In these studies, the reconstructions of the amount of SiH2 relied on the assumption of the Boltzmann thermal distribution when estimating the population of the lower state. No account of the possible non-LTE population of SiH2 molecules after dissociation was made, which could potentially hamper the count of the SiH2 molecules and affect the reaction rates estimated. A similar experimental technique was used in Hertl and Jolly71 to monitor SiH2 in SiH4 plasma.</p><p>It is the second, dissociation, part of the reaction shown in Fig. 1 (Si2H6 → SiH2 + SiH4) we study in this work. More specifically, we show that (i) the (vibrational) populations of the molecules produced in reactions can be very different from the Boltzmann distribution and is important to take into account when interpreting spectroscopic measurements. That (ii) spectral shapes of the dissociated SiH2 can bear strong non-LTE character, very different from the LTE spectrum of an LTE SiH2 sample making it possible to distinguish between different reaction stages and even between different dissociation channels the silylene molecules it is produced from. In this work non-LTE spectra of SiH2 under conditions similar to dissociation processes expected in these experiments are modelled.</p><p>Recently we have computed an accurate ro-vibrational line list for SiH2, named CATS.72 It covers a large range of rotational and vibrational excitation, capable of modelling very hot spectra of this molecule (up to T = 2000 K) as part of the ExoMol database.73 The CATS line list was produced using the program TROVE,74,75 which solves the nuclear motion Schrödinger equation variationally. The ro-vibrational energies and corresponding wavefunctions were computed using an accurate, empirically refined potential energy surface (PES) of silane and a high-level ab initio dipole moment surface (DMS). The ro-vibrational probabilities (in the form of Einstein A coefficients) were computed using a high level ab initio dipole moment surface.</p><p>The study by Clark et al.72 forms the basis for the present work, where we utilize the CATS line list, wavefunctions, purpose-built numerical basis set, and the CATS computational TROVE setup to model non-LTE spectroscopic properties of SiH2 produced from dissociation of Si2H6 through different reaction channels. Using a simplified 1D Harmonic oscillator wavefunctions (see Pastorek et al.6) and more sophisticated 3D vibrational CATS wavefunctions from accurate variational calculations, the non-LTE ro-vibrational populations of SiH2 are generated and used to produce non-LTE spectroscopic spectra of different dissociation channels of disilane. To this end we investigate reaction topology connecting the global minimum of Si2H6 with the closest saddle points and local minima as well as the corresponding structural properties using a high level ab initio theory cc-pVTZ-F12/CCSD(T)-F12b76,77 employing the program MOLPRO2015.78</p><p>Theoretically, the non-LTE properties of dissociating molecules were studied by Band and Freed.79 In the present work we use general approach of Berry80 and Band and Freed,79 which assume no significant structural changes between the reactant and product nuclear configuration, along with the slow vibrational relaxation of the product13 to investigate non-LTE spectroscopic signatures of SiH2 produced from dissociation of Si2H6. In this paper we specifically consider situations where the vibrational relaxations are not achieved during the time of the experiment, so that the molecules still hold the memory of the structure during the reaction or dissociation. The rotation relaxation time however is much shorter and the rotational degrees of freedom can be usually assumed to satisfy the Bolzmann equilibrium.13</p><p>We also investigate possible non-LTE impact on the electronic Ã(0,2,0)–X̃(0,0,0) spectrum of SiH2. This is a favorite spectroscopic system for the detection of SiH2 due to the large Franck–Condon factor and the availability of suitable laser.39,70,81–91</p><p>The non-LTE absorption spectra of SiH2 are simulated using the (non-LTE) ExoCross program,92 where a new feature of non-Boltzmann populations was added. ExoCross has been previously used to model spectra of molecules in environments that can be characterized using two temperatures, vibrational and rotational.93–96</p><p>The paper is structured as follows. In Section 2 we describe the calculations of potential energy surfaces for the disilane and silylene structures. The theory used in this paper is described in Section 3. In Section 4.1 we calculate the 1D harmonic wavefunction population and use them to produce non-LTE spectra of SiH2 corresponding to different dissociation routes. In Section 4.2 we calculate the populations and subsequent non-LTE ro-vibrational spectra of SiH2 using the full 3D wavefunctions and describe the new TROVE methodology. A non-LTE electronic Ã(0,2,0)–X̃(0,0,0) spectrum of SiH2 is presented in Section 4.4. Conclusions are offered in Section 5.</p><!><p>In order to better understand the reaction process of breaking Si2H6, the topology of Si2H6 has been investigated by performing a structural analysis of Si2H6 using a high level ab initio theory. This includes finding the global minimum (GM), local minima (LM), transition states (TS), reaction barriers as well as reaction paths, as detailed below. A reaction slice through the global PES of Si2H6 helps to indicate how likely local minima or transitions states were to be formed based on corresponding topology. These properties of Si2H6 were obtained using the geometry optimization and reaction path finder implemented in MOLPRO201578 using the explicitly correlated coupled cluster method CCSD(T)-F12b76,77 with the F12-optimized correlation consistent basis set, VTZ-F1297 in the frozen core approximation. The calculations employed the diagonal fixed amplitude ansatz 3C(FIX)98 and a Slater geminal exponent value of β = 1.0a0−1.99 The auxiliary basis sets were chosen to be the resolution of the identity OptRI100 basis and the aug-cc-pV5Z/JKFIT101 and cc-awCV5Z/MP2FIT102 basis sets for density fitting. In the following this level of theory will be referenced to as VTZ/CCSD(T)-F12b.</p><p>We shall refer to different disilane isomers as dGM, dLM and dTS to distinguish them from the SiH2 fragments GM, LM and TS as discussed below.</p><p>The global minimum of SiH2 (dGM) has a symmetrical, staggered D3d structure. The closest local minimum (dLM) has an inverted, C3v structure. The lowest transition state, which will be referred to as dTS (TS1 in Becerra et al.46 and TS2 in Tonokura et al.103), has also inverted structure, just a few kJ mol−1 above dLM.46 These structures together with the corresponding optimized parameters for three geometries most relevant to our work are collected in Table 1. Our structural parameters compare well with that from the literature. The structure of the dGM has also been determined spectroscopically, with the equilibrium bond lengths rSi–H = 1.492 Å, rSi–Si = 2.331 Å and bond angles β∠HSiSi = 110.3° and α∠HSiSi = 108.6°.104 The structure of the deuterated isotopologue Si2H5D was reported as rSi–H = 1.4874(17) Å, rSi–Si = 2.3317(15) Å, β∠HSiSi = 110.66(16)°.105</p><p>The reaction path connecting the disilane isomers dGM, dLM and dTS is shown in Fig. 2. A zoom of the dLM side is shown as inset. The energies of the global and local minima are 178 kJ mol−1 and 4.1 kJ mol−1 below the transition state, respectively. The energy and geometry information for dLM, dGM, dTS are collected in Table 2. The results calculated compare well with both the results of Becerra et al.46 and Sakai and Nakamura,60 albeit both vary for the dGM structure by 36 kJ mol−1.</p><!><p>TS2 in Sakai and Nakamura60 is the same as TS1 in Becerra et al.46</p><p>Complex 2 in Sakai and Nakamura60 is the same as LM1 in Becerra et al.46</p><p>Relative to SiH2 + SiH4, kJ mol−1.</p><p>Relative to dTS.</p><!><p>The dLM isomer of Si2H6 has a shallow potential with a very low barrier to dTS of 344 cm−1. It can be also recovered using lower levels of theory, for example, using MP2/6-311G(d,p)46 and even with the UFF force fields implemented in Avogadro 1.2.0106via the steepest descent method and 4 steps per update.</p><!><p>With the aim to give more quantitative information on the structural and dynamical properties of five SiH2 fragments from Si2H6, the fragments are described as follows. GM is an SiH2 fragment from the global minimum structure (dGM); LM-L is an SiH2 fragment from the left hand side (LHS) of the local minimum structure (dLM) with the Si–H3 umbrella group pointing outside; LM-R is an SiH2 fragment from the right hand side (RHS) of the dLM structure with Si–H3 pointing inside; TS-L is an SiH2 fragment from the LHS of the transition state structure (dTS), Si–H3 umbrella group points outside; and TS-R is an SiH2 fragment from the RHS of the dTS structure, Si–H3 umbrella group pointing inside. The structural parameters and structures are shown in Table 3 and the energies are shown in Table 2. The harmonic frequencies for the disilane molecules were computed with MOLPRO using the TVZ/CCSD(T)-F12b level of theory. The results from the frequency analysis are shown in Table 3 along with a comparison with literature (both experimental, if available, or theoretical data is shown). The symmetry of dGM is D3d whereas the symmetry for the dLM and dTS are both C3v. The columns titled "Theory" were calculated in this work, and the degenerate states have been removed, with an average calculated if there were any differences owing to computational errors associated with the lower symmetry used by MOLPRO. All frequencies of dLM are positive thus confirming that it is a minimum with a stable structure. The 'negative' (or imaginary) harmonic frequency of dTS is −551.6 cm−1.</p><!><p>Confirmed by the experimental gas phase Raman spectrum of Durig and Church.59</p><p>Although there are three lines of 925 cm−1, there are no triply degenerate states in the C3v symmetry. Roenigk et al.67 assign one 925 cm−1 to ν15, which in the dGM symmetry is an E state degenerate with ν16. The dTS ν16 can be found at 960 cm−1. The other two 925 cm−1 line were assigned to ν9 and ν10 (E).</p><!><p>We now consider the decomposition reaction Si2H6 → SiH2 + SiH4 shown in Fig. 1 and model the vibrational populations of the product SiH2 assuming that the corresponding relaxation time to LTE is much longer than the time of the spectroscopic experiment. We aim at simulating non-LTE IR spectra of SiH2 using these populations to demonstrate their impact on the spectral shape of dissociated species.</p><p>In our description of the non-LTE vibrational population of the dissociated molecule we follow the polyatomic Franck–Condon type approximation by Band and Freed79 and Berry80 based on the structural differences between reactant and product assuming no significant change in nuclear configuration of the molecule. In order to connect the product (gas phase SiH2) to an initial structure of Si2H6 through the dissociation process, we assume that the dissociation happens instantaneously, i.e. the initial configuration of the product SiH2 corresponds to the structural parameters (bond lengths Si–H and inter-bond angles ∠HSiH) of SiH2 as a fragment of Si2H6, for which the parameters collected in Table 1 are used. For example, for the dissociation from the dGM structure, the initial configuration of the gas phase SiH2 is assumed to be rSi–H = 1.482 Å, α∠HSiH = 108.6°. Naturally, this is a very deformed geometry comparing to the equilibrium structure of the gas phase SiH2, re = 1.5137 ± 0.003 Å and αe = 92.04 ± 0.05°.107 After being dissociated, in relaxing to be a free molecule, the fragment SiH2 has added vibrational energy and is hence in a situation when it's vibrational populations do not match the LTE distribution for corresponding temperature of the surroundings, at least for the vibrational degrees of freedom.</p><p>Experience shows that the rotational equillibration time is usually very short and we can thus safely assume the LTE conditions for the rotational degrees of freedom with the rotational temperature the same as the temperature of the surroundings (see also e.g. Dudás et al.13). The vibrational population however is not in the thermal equilibrium and therefore no sensible vibrational temperature could be associated with the corresponding population.</p><p>The non-LTE population NJ,k,ν(T) of a ro-vibrational state |J, k, ν〉 is then given by:61where Nvibν is a non-LTE vibrational population, Trot is the rotational temperature, J is the total angular momentum quantum number; k is a generic rotational quantum number, e.g. the projection of the total angular momentum on the molecular z axis; is a generic vibrational quantum number/label, e.g. a combination (ν1,ν2,ν3) to describe vibrational states of a triatomic molecule; gnsi is the nuclear spin degeneracy; T is the temperature, c2 is the second radiation constant. In eqn (1), ẼrotJ,k is the rotational part of the ro-vibrational energy term value approximated as2ẼJ,k,ν = Ẽvibν + ẼrotJ,k,where hc·Ẽvibν is the vibrational (J = 0) energy ('band center') and hc·ẼJ,k,ν is the total ro-vibrational energy. The non-LTE partition function QnLTE in eqn (1) is given by3Our aim is to calculate the vibrational populations Nv ≡ Nvibν of SiH2 as produced by instantaneous (vertical) dissociation from three structures, dGM, dLM and dTS. In case of dLM and dTS, the SiH2 fragment can originate from any of the two different sides of their inverted structures, which should be taken into account. We therefore have to consider five different fragments, as shown in Table 4.</p><p>Let us assume that the Si2H6 is LTE and hence is in its ground vibrational state at the moment of dissociation, while the fragment SiH can end up in any vibrationally excited state ≡ (ν1,ν3) with some transition probability giving rise to the vibrational population Nvibν. On top of that we also assume a full separation of the stretching Si–H and bending ∠HSiH modes inside Si2H6 in its ground vibrational state. Possible consequences of deviation from these approximations are discussed below.</p><p>Under the assumptions made we define the vibrational population of the SiH2 fragment as a Franck–Condon factor for a vertical transition from the ground vibrational state of disilane Si2H6 to SiH2 + SiH4 with the gas phase (g.ph.) SiH2 transferred to some vibrational state |(g.ph.)〉 = |ν1,ν2,ν3〉. In the approximation of the full separation of the fragment SiH2 from the rest of Si2H6, the population of SiH2 can be represented as an overlap between the ground state wavefunction |ν′′ = 0(fragment)〉 of a fragment SiH2 and that of the corresponding vibrational state |ν′ = (g.ph.)〉 of the gas phase SiH2 as given by4N = |〈0(fragment)|ν(g.ph.)〉|2.Here the spacial wavefunction of the SiH2 fragment (i.e. a combination of two adjacent Si–H bonds in disilane with an angle between them forming the dissociating SiH2) is projected on vibrational eigenfunctions of the gas phase SiH2 to give the corresponding populations of SiH2. The highest populated energy level will have the largest overlap between these wavefunctions.</p><p>The calculated temperature dependent populations NJ,k,ν(T) in eqn (1) can be then combined with a molecular line list for SiH2 to simulate absorption or emission spectra of this molecule under the non-LTE conditions as defined by NJ,k,ν(T) in eqn (1). Here we use the ExoMol line list CATS by Clark et al.72 as provided by ExoMol (www.exomol.com). Technically this is done by incorporating the non-LTE vibrational densities N into the ExoMol States file as described in Section 4 (the ExoMol file formats are discussed extensively elsewhere73). A non-LTE spectrum of SiH2 for given T and P is then calculated using CATS' Einstein A coefficients with the ExoCross program,92 where a new non-LTE option has been implemented as part of this work. The rotational populations are assumed to be in LTE according with eqn (1).</p><!><p>Two approaches were used for the calculation of the population densities of the fragment SiH2. One approach – named the decoupled 1D approach – is where the 3D wavefunctions of the fragment as well as of free SiH2 are represented by products of 1D parts with the harmonic oscillators as wavefunctions. This simplified model is mainly used to illustrate the idea of our non-LTE treatment. The second, more accurate approach – named the 3D approach – is based on the full 3D vibrational wavefunctions computed using the variational program TROVE.74 Both approaches are presented in the following in order to assess and compare the accuracy achieved.</p><!><p>A vibrational state |ν〉 = |ν1,ν2,ν3〉 of SiH2 is characterized by the three (normal mode) quantum numbers ν1, ν2 and ν3 corresponding to the two stretching modes (ν1 and ν3) and one bending mode (ν2) of the SiH2 molecule. The 1D approach considers the stretching Si–H1, Si–H2 and bending ∠HSiH modes, both of the molecular fragment and gas phase SiH2 molecules, as fully independent and described by one dimensional (1D) wavefunctions under the harmonic approximation, as given by:5Here x is a dimensionless coordinate describing either the stretching r = rSi–H or bending α = α∠H–Si–H coordinate as follows:67with89and ωs = ωSi–H, ωb = ω∠HSiH. In eqn (5)H(x) is a Hermite polynomial and Cv is the corresponding normalization constant.</p><p>The constants as and ab correspond to inverse masses of the vibrational part of a free three-atomic molecule expressed in terms of the internal coordinates r1, r2, α (see, e.g. Sutcliffe and Tennyson,108 Yurchenko et al.109).</p><p>A 1D population for the corresponding mode of the gas phase SiH2 molecule is given by eqn (4) with |ν〉 = Ψv(x). The different disilane fragments have different structural parameters re, αe, ωs and ωb, see Tables 1 and 3, and thus lead to different ground state vibrational 1D wavefunctions |0(fragment)〉 (stretching or bending) and hence result in different vibrational populations Nν of the gas phase SiH2 according with eqn (4).</p><p>For each of the three modes (two stretching and one bending), 1D wavefunctions of the gas phase SiH2 for 30 vibrational states from ν0 up to ν29 were calculated. These 1D wavefunctions were then numerically integrated with the corresponding ground state 1D wavefunctions |ν = 0〉 of the fragment in question.</p><p>The total vibrational population Nν1,ν2,ν3 in this approximation is then given by a product10</p><p>where  and  are obtained using the stretching harmonic oscillators wavefunctions |ν1〉 and |ν3〉, while  is obtained using the corresponding bending harmonic oscillator wavefunction |ν2〉. The independent treatment of the two stretching populations is partly justified by the local mode character of the vibrational degrees of freedom of SiH2 due to the near 90° bond angle (see, e.g. Jensen110 and Clark et al.72). The asymmetric vibrational modes of SiH2 (B2 in C2v) are non populated in this 1D approximation. This is because for the parallel nature of the Franck–Condon transitions from the ground vibrational state of disilane, which is fully symmetric (A1), only symmetric states of SiH2 give rise to non-zero integrals in eqn (4). For example, the vibrational population of the (2,1,0) state (2 stretching and 1 bending quanta of A1) of the GM fragment is obtained as a product of ,  and , respectively, resulting in N(ν2,ν1,ν0) = 0.088, while the population of the B2-type (0,0,1) vibratitonal state is assumed to be zero. The populations Nν1,ν2,ν3 are pre-calculated for each vibrational state (ν1,ν2,ν3) of SiH2 and added to the CATS State file to be used in non-LTE simulations (see below, Section 4.3). All the vibrational populations computed and used as part of this work are provided in the ESI.†</p><!><p>Examples of overlapping bending mode Harmonic wavefunctions used in calculations of populations for ν = 0, 1 and 2 of the five fragments are shown in Fig. 3 and 4 for the bending and stretching modes respectively. In all cases the black curve represents the ground state wavefunction for the non-LTE fragment, and the blue, green and red line show the ν = 0, ν = 1 and ν = 2 LTE wavefunctions of gas phase SiH2.72 The corresponding 1D populations as an integral of the overlaps between the LTE and non-LTE wavefunctions are plotted in Fig. 5 for the bending and stretching modes of the five fragments. The stretching populations exhibit a typical Boltzmann-like distribution with the ground vibrational state ν = 0 as the mostly populated in all five cases. This is expected because their equilibrium bond lengths are rather similar to that of the gas phase SiH2. In case of the bending populations, only LM-R and TS-R have ν = 0 to be with the highest populations, while for LM-L, TS-L and GM the N distributions exhibit strong non-LTE character with ν = 2 to be almost as populated as the ground vibrational state ν = 0.</p><p>The shapes and positions of the curves in Fig. 3 and 4 match with the parameters from Table 1. The larger ω used for the bending modes manifests itself as wider curves, while the black curves are all centred around the equilibrium bond angles and length listed in Table 1.</p><p>The similarity between the populations of the LM-L, TS-L and GM fragments is expected owing to their similar structural parameters. It is interesting to see the most populated vibrational levels of LM-R and TS-R are always lower than the vibrational levels of LM-L, TS-L and GM.</p><p>Fig. 6 and 7 show the 1000 cm−1 and 2000 cm−1 (10 μm and 5 μm) bands for the SiH2 absorption spectrum, respectively, simulated using the non-LTE densities from Fig. 5 for all five cases considered and compared to the LTE scenario assuming the (rotational) temperature of T = 296 K and using the CATS line list. The strongest bands are indicated using different colours. Fig. 6 focuses on the 1000 cm−1 band. Most of the non-LTE spectra contain bending hot bands (020)–(010) and (020)–(010), which are stronger than the fundamental band (000)–(000). It can be seen that the P and R branches of the non-LTE spectra are shifted to lower wavenumbers in the GM, LM-L and TS-L spectra. In the TS-R and LM-R spectra the bands are not shifted, with the LM-R spectrum having only the fundamental (010)–(000) band visible.</p><p>The plots in Fig. 7 show the 2000 cm−1 band in the region of the polyad (100)/(020)/(001) for the five fragments, with the strongest fundamental band (001)–(000). The non-LTE intensities of the hot bands (011)–(010), (200)–(020) are found to be comparable to the intensities of the (001)–(000) band. The Q branch is clearly shifted for the GM, TS-L and LM-L molecules. The band is less shifted for the TS-R and LM-R fragments, but owing to the increased similarity between the fragment and molecular structures with TS-R and LM-R this is to be expected. Only the main polyad system (100)/(020)/(001) is visible for the TS-R spectrum (indicated as (001)–(000) in Fig. 7).</p><p>With the equilibrium structures of the TS-R and LM-R fragments being similar to the equilibrium structure of SiH2, their non-LTE spectra are expected to be a similar spectrum to LTE. Indeed, for the 1D harmonic approach their P, Q and R branches maintain the expected LTE intensities for both the 1000 cm−1 and 2000 cm−1 bands.</p><!><p>In a full 3D approach, the ground state wavefunction |0,0,0(fragment)〉 represents an SiH2 fragment of an 18D ground state vibrational wavefunction of Si2H6:ϕSi2H6 = |0,0,0(SiH2)〉 |0,0,0,0,0,0,0,0,0,0,0,0,0,0,0(SiH4)〉.Here we assume the approximation that the corresponding three modes (Si–H1, Si–H2 and ∠H1SiH2) are independent from the rest of the molecule so that all other modes, not relevant for the gas phase SiH2, can be eliminated (integrated out), including the reaction coordinate and vibrational modes of SiH4. This is in line with the assumptions used previously by Band and Freed79 and Berry.80 Apart from this approximation we will treat the SiH2 fragment as accurate as possible. The corresponding wavefunction |0,0,0(fragment)〉 is obtained by solving a 3D vibrational Schrödinger equation for these three degrees of freedom with a realistic PES obtained using a high level of ab initio theory (the same as above, VTZ/CCSD(T)-F12b with MOLPRO).</p><p>The vibrational populations are then modelled using the Franck–Condon integrals as follows:11Nν1,ν2,ν3 = |〈0,0,0(fragment)|ν1, ν2, ν3(g.ph.)〉|2,with the sum of all populations over all states totalling 1. In this equation, |ν1,ν2,ν3(g.ph.)〉 is an accurate vibrational wavefunction of a gas phase SiH2 molecule, obtained by solving the vibrational Schrödinger equation with an accurate PES. We use the TROVE variational program and the refined PES of SiH2 by Clark et al.72 to generate |ν1,ν2,ν3(g.ph.)〉 for all vibrational excitations required. For the details on the TROVE calculations see below and also Clark et al.72</p><p>In order to simplify the 3D integration in eqn (11), the variational wavefunction |0,0,0(fragment)〉 is obtained using the same vibrational basis set as the variational solution of the gas phase SiH2. By taking advantage of the compatibility of the orthogonality of the basis sets, the Franck–Condon factors are then given by as a sum of products12of the corresponding eigen-coefficients Ci1,i2,i3(fragment) and Ci1,i2,i3(g.ph.), obtained variationally in independent calculations using a new implementation in TROVE.</p><p>TROVE uses optimized non-standard vibrational basis sets, generated numerically by solving 1D Schrödinger equations for realistic 1D potentials.74 This procedure allows producing compact basis functions optimized for a specific problem. In our case, the PESs of the corresponding five fragments and of the gas phase SiH2 are different and therefore the generated basis sets would be different and even not orthogonal. We therefore implemented a feature in TROVE allowing to read and use externally generated basis functions. Of course all relevant calculation setups must be compatible, including the numerical grids used for the stretching and bending modes and their sizes. Using foreign basis sets certainly degrades their quality. However, since we are only interested in fragments' ground state wavefunctions, this degradation can be mitigated by including enough basis functions. Our typical 1D basis sets contain 12–24 functions (see details below), which should be more than enough to obtain a converged ground state solution even with non-optimized basis sets.</p><!><p>For our new 3D populations corresponding to dissociations from the fragments, five PESs were generated as follows. We assume that PESs of a dissociating Si2H6 molecule can be approximated as a sum of two independent fragments:VSi2H6 = VSiH2(r1, r2, α) + Vrest,where the individual stretching and bending modes of SiH2 fragments are fully separable:13VSiH2(r1, r2, α) = fstr(r1) + fstr(r2) + fbnd(α).The stretching part of the potential is given by a Morse-like expansion14while the bending part is a Taylor-type expansion in terms of the displacement from the corresponding equilibrium value:15The expansion constants fi and gi representing fragments' potential energies VSiH2(r1,r2,α) were obtained by fitting eqn (14) and (15) to the ab initio data computed as 1D slices on the global surfaces for the five fragments from Si2H6 (dGM, dLM and dTS) using VTZ/CCSD(T)-F12b consisting of 24 bending and 34 stretching geometries, distributed around the corresponding equilibria. They are listed in Table 5.</p><p>Fig. 8 illustrates the ab initio PESs of different fragments as 1D cuts for the stretching and bending modes compared to the corresponding cuts of the gas phase SiH2 molecule. The bending cuts have especially different shapes with shifts to larger equilibrium angles and much steeper PESs. The differences in the stretching cuts are less pronounced. These features are important for the non-LTE behaviour of the corresponding excited states populations, with the bending degree of freedom to have stronger non-LTE character than stretching.</p><!><p>The vibrational wavefunctions of SiH2 were computed using the variational nuclear motion program TROVE with the same setup as in Clark et al.72 Details of the TROVE methodology are discussed extensively elsewhere.74,111–113 Here, we give a brief outline of the main calculation steps. The TROVE kinetic energy operator is Taylor expanded up to sixth order around the SiH2 equilibrium geometry in terms of linearized coordinates.114 The primitive basis set is constructed from 1D mode numerical basis functions using the Numerov–Cooley approach115,116 by solving three 1D Schrödinger equations, for each vibrational degree of freedom. The stretching basis functions are then improved by solving a 2D Schrödinger equation for a reduced stretching Hamiltonian. The resulting stretching eigenfunctions are contracted, classified according with the C2v(M) symmetry group114 using an optimized symmetrization procedure112 and combined with the bending primitive basis functions to form our final, symmetry-adapted 3D vibrational basis set. The basis set coverage is defined by the polyad number cut-off16P= 2(ν1 + ν3) + ν2 ≤ 24,where ν1 and ν3 are the stretching and ν2 is the bending quantum numbers, with the maximal excitations 12, 12 and 24 respectively.</p><p>For the gas phase SiH2 calculations we employed the empirically refined PES by Clark et al.72 For the five SiH2 fragments only the vibrational ground state wavefunctions |0,0,0〉 were computed using the same setup and utilizing the basis functions from the g.ph. calculations as described above, but for the fragments' ab initio PESs from eqn (13).</p><!><p>The 3D vibrational populations Nν1,ν2,ν3 in eqn (11) were computed using the corresponding TROVE eigenfunctions viaeqn (12).</p><p>The 3D results generally agree with the 1D case. As an example, Fig. 9 shows the TS-L case: the bending primitive wavefunctions of the gas phase SiH2 are compared to the ν = 0 bending wavefunction generated for the PES of TS-L. Fig. 10 shows the corresponding non-LTE vibrational populations as a function of the corresponding energies for the same TS-L scenario of the 3D vibrational populations. These populations are very different from the Boltzmann distribution, also shown on this figure for the SiH2 vibrational states at T = 1000 K, which exhibits an exponential decay with ν = 0 at its maximum. It nicely demonstrates that it would not be possible to associate a single vibrational temperature for this bell-shaped distribution.</p><!><p>A non-LTE absorption line intensity Ifi (cm per molecule) can be calculated as17where Afi is the Einstein A coefficient (s−1), 


<svg xmlns="http://www.w3.org/2000/svg" version="1.0" width="13.454545pt" height="16.000000pt" viewBox="0 0 13.454545 16.000000" preserveAspectRatio="xMidYMid meet"><metadata>
Created by potrace 1.16, written by Peter Selinger 2001-2019
</metadata><g transform="translate(1.000000,15.000000) scale(0.015909,-0.015909)" fill="currentColor" stroke="none"><path d="M160 840 l0 -40 -40 0 -40 0 0 -40 0 -40 40 0 40 0 0 40 0 40 80 0 80 0 0 -40 0 -40 80 0 80 0 0 40 0 40 40 0 40 0 0 40 0 40 -40 0 -40 0 0 -40 0 -40 -80 0 -80 0 0 40 0 40 -80 0 -80 0 0 -40z M80 520 l0 -40 40 0 40 0 0 -40 0 -40 40 0 40 0 0 -200 0 -200 80 0 80 0 0 40 0 40 40 0 40 0 0 40 0 40 40 0 40 0 0 80 0 80 40 0 40 0 0 80 0 80 -40 0 -40 0 0 40 0 40 -40 0 -40 0 0 -80 0 -80 40 0 40 0 0 -40 0 -40 -40 0 -40 0 0 -40 0 -40 -40 0 -40 0 0 -80 0 -80 -40 0 -40 0 0 200 0 200 -40 0 -40 0 0 40 0 40 -80 0 -80 0 0 -40z"/></g></svg>

fi is the transition wavenumber (cm−1), QnLTE(T) is the non-LTE partition function defined in eqn (3).</p><p>In order to simulate absorption spectra of the gas phase SiH2 assuming a non-LTE vibrational populations Nν1,ν2,ν3, we use the line positions and Einstein A coefficients from the ExoMol CATS line list72 employing the ExoCross program.92 The ExoMol line lists are formatted as two files, a States file and a Transition file. It is described extensively elsewhere73 and in this paper we shall only discuss how the ExoMol format has been adapted for use in non-LTE situations.</p><p>To adapt the CATS States file72 for non-LTE applications an additional 'density' column Nν1,ν2,ν3 was added as a final column, see an extract from the States file in Table 6. This column contains the weightings to the transition probabilities as populations of the vibrational levels occupied by the gas phase SiH2. The 'density' column is specific for the calculation of line lists for non-LTE molecules is not routinely included into the ExoMol States files.</p><p>This column is read by ExoCross and used to give the population weighting to each line intensity. Evaluation of the non-LTE population NJ,k,ν(T) given for each ro-vibrational state as in eqn (1) is based on the knowledge of the corresponding vibrational state (ν1,ν2,ν3) as well as the rotational energy contribution ẼrotJ,k. Therefore for this approach to work it is mandatory for all ro-vibrational states to be vibrationally assigned in order to be able to subtract the vibrational contribution Ẽvvib from the total energy according with eqn (2). All vibrational quantum numbers (ν1,ν2,ν3) are not required, only a vibrational index indicating the vibrational state in question. In our model, all asymmetric vibrational states (B2) are not populated due to the zero overlap with the ground state of the A1 symmetry in eqn (11), as part of the completely vertical Franck–Condon approximation.</p><!><p>Fig. 11 shows non-LTE spectra of SiH2 in the two main spectroscopic regions, 1000 and 2000 cm−1 for all five fragments considered (GM, TS-L, TS-R, LM-L and LM-R). The corresponding vibrational non-LTE populations were generated with the 3D TROVE approach. The rotational populations assume the Boltzmann distribution with the rotational temperature of Trot = 296 K. These non-LTE spectra are compared to the corresponding LTE spectra of SiH2 at T = 296 K, which comprise mainly two fundamental bands, (0,1,0)←(0,0,0) (1000 cm−1 region) and (0,0,1)←(0,0,0) (2000 cm−1 region), with the hot bands suppressed due to the relatively low temperature. The non-LTE spectra are dominated by the hot bands (0,2,0)←(0,1,0), (0,3,0)←(0,2,0) (1000 cm−1 band), (1,1,0)←(0,1,0) and (1,1,0)←(0,2,0) (2000 cm−1 band). The centers of the hot bands are systematically red shifted compared to the fundamental band centres and serve as distinct signatures of the non-LTE effects. The Q-branch of (011) ← (020) is especially distinct compared to the Q-branch of the LTE fundamental band in the 2000 cm−1 region, with the difference of about 15 cm−1. Only the LM-R spectra of non-LTE are very similar to the LTE spectra. This is not surprising considering that the equilibrium values of rSiH and α∠H–Si–H of the LM-R structure are very similar to the corresponding equilibrium parameters of the gas phase SiH2.</p><p>The stark differences of the non-LTE spectra offer an ability for experiment to distinguish between SiH2 molecules produced from fragmenting Si2H6, and even to indicate dissociation channels involved.</p><!><p>The visible electronic band Ã–X̃ of SiH2 has often been used to study different reactions involving leading to silylene.91 Here we used the program RENNER117 to simulate a non-LTE electronic spectrum of SiH2 with the spectroscopic model by Yurchenko et al.118 for the Ã1B1–X̃1A1 system. The model includes two empirically adjusted PESs, for the Ã and X̃ states, respectively, and an ab initio (MRCI) Ã–X̃ transition dipole moment surface (TDMS).</p><p>In the RENNER calculations, the size of the basis set originally used in Yurchenko et al.118 was reduced in order to be able to increase the rotational excitations. The main purpose of this exercise is to show a qualitative impact of the non-LTE populations on the spectral shape of the electronic band and not so much the quality of the line positions, and therefore a smaller basis set is justified. We used 18 and 12 bending basis functions for the X̃ and Ã electronic states, respectively, for every |k| block (where |k| ≤ J). The X̃1A1 electronic state basis set included NA = 12 stretching functions of the A1 symmetry and NB = 10 stretching functions of the B2 symmetry. For the Ã1B1 state NA = 10 and NB = 8 stretching functions were used. These stretching functions were constructed from the Morse oscillator functions |n1〉|n3〉 with n1 + n3 ≤ Nstretch = 12.</p><p>A rovibronic line list for the Ã1B1–X̃1A1 of SiH2 was generated covering the rotational excitations up to Jmax = 15 with the lower state energies (X̃) truncated at hc·25 000 cm−1 and the upper state energies truncated at hc·28 000 cm−1.</p><p>For the non-LTE simulations we used the 1D vibrational population model with the structural parameters corresponding to TS-L from Table 3. Fig. 12 shows a non-LTE electronic spectrum of SiH2 in the region of the band Ã(0,2,0)←X̃(0,0,0), assuming the rotational temperature Trot = 500 K, compared to an LTE spectrum of T = 500 K. The non-LTE spectrum contains the hot band Ã(0,3,0)←X̃(0,1,0) which can be used to identify the non-LTE character of the system. The rovibronic line 101←110 belonging to this band used in a number of experimental studies involving SiH2 as a reaction product71,89,90 to estimate reaction rates. It is common for such studies to assume the Boltzmann equilibrium at different stages of the analysis of the measurements. For example, the partition function of SiH2 is required to estimate the number density of SiH2 in its lower, ground electronic state,70,71,119 which is directly affected by the LTE assumption. As our calculations show, the number densities of SiH2 as a reaction product can vary significantly under non-LTE depending on the reaction pathway and impact the experimental rates.</p><!><p>The focus of this paper is on the new features which have been added to TROVE to allow modelling the non-LTE populations of polyatomic molecules. We have demonstrated this capability by modelling the non-LTE line list of SiH2 calculated with 3D wavefunctions and TROVE, and compared them to non-LTE spectra of SiH2 modelled using a 1D harmonic approach, and the LTE line list calculated previously by ExoMol.</p><p>There are two stable isomers of disilane, a local minimum structure and a global minimum structure, with a third transition state structure also known. Non-LTE spectra of SiH2 corresponding to dissociation of disilane from different sides of the three disilane isomer were computed. We have shown that the non-LTE spectra of SiH2 are different in most cases. This is important as the spectrum of SiH2 is used to monitor the quantity present in a reaction as a means to track the progress of SiH2 + SiH4 → Si2H6 and Si2H6 → SiH2 + SiH4 when calculating the corresponding rate constant. If Si2H6 is decomposing at a rate slower than it is being formed, then tracking the quantity of SiH2 can give a rate constant that is not reflective of the speed of reaction, and merely an indication of the equilibrium balance of the two species SiH2 and Si2H6.</p><p>In two approaches considered, 1D and the 3D, we assume that the rotational degrees of freedom are equillibrated quickly once the dissociation from disilane occurs, hence we use the Boltzmann distribution for the rotational degrees of freedom. We also assume that the SiH2 fragment during the instantaneous dissociation and is fully decoupled from the rest of the Si2H6 molecule, i.e. can be described by a 3D wavefunction in its lowest, relaxed vibrational configuration and has the same structural parameters as the Si2H6 molecule.</p><p>We have shown that the non-LTE spectra of SiH2 can be calculated by the new TROVE methodology and existing ExoMol line list, and it compares well to the simpler 1D harmonic approximation published previously. The method could be applied to the non-LTE spectroscopy of other small molecules including SiH4, which has not been explored here.</p><p>We have also shown that despite the many approximations used in the 1D approximation (separability of the modes, Harmonic approximate etc.), the results compare well to the results obtained using the full 3D approach. This lends confidence in using the simplified but robust 1D approach in similar non-LTE studies, as e.g. we have used to model the CO non-LTE spectra,6 recently, which are planning to explore in the future.</p><p>The methods described here can be used model the intensity distribution of the reaction products and to ascertain from what molecule the SiH2 dissociated from. The equilibrium structured parameters (bond lengths and angles) can be treated as effective parameters to be adjusted to reproduce the experimental spectra.</p><!><p>There are no conflicts to declare.</p>