Modelling the non-local thermodynamic equilibrium spectra of silylene (SiH₂)

Abstract

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 (SiH₂) produced in a decomposition of disilane (Si₂H₆), a reaction of technological importance. Two approaches for non-LTE vibrational populations of the product SiH₂ are introduced: a simplistic 1D approach based on the Harmonic approximation and a full 3D model incorporating accurate vibrational wavefunctions of SiH₂ 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.

---

1 Introduction

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

The modern day study of non-LTE spectroscopy can be traced as far back as the 1930s to the original papers of Milne,¹⁴ 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.,²⁵ Funke et al.,²⁶ Pannier and Laux²⁷ and references therein. As such, non-LTE spectra are often vital for the modelling of astrophysical problems, including planetary atmospheric properties,²⁸ stellar atmospheres of solar system and exoplanets^29–36 and the ISM.^25,37,38

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.² 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.² 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 (SiH₂) produced from disilane (Si₂H₆).

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 Si₂H₆, are just as important to study as the elementary reactions.⁴² Silane containing reactions are also of importance for astrophysics, with the presence of SiH₄ in IRC + 10216 discussed by Goldhaber and Betz,⁴³ Kaiser and Osamura.⁴⁴

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 Si₂H₆e.g.,^46–53 as well as of the formation and reactions involving this molecule.^54–58 Agrawal et al.⁵³ 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.⁵⁰ reported a force field for Si₂H₆. The main isomer of Si₂H₆ has a staggered, ethane-like, structure (see Fig. 1) with a low barrier (1.2 kcal mol⁻¹) to the eclipsed conformation.⁵⁹ Si₂H₆ has been shown to have a local minimum as an inverted stable structure with one of the Si–H₃ ‘umbrellas’ pointing to the center as well as a transition state with a similar inverted configuration.^46,60 These structures are nominally asymmetric (C_s symmetry) but essentially acquiring the C_3v symmetry.

Fig. 1 Schematic reaction path for the silylene. Molecular SiH₂ is the starting material, fragment SiH₂ is the product.

Thermal decomposition of Si₂H₆ has been extensively studied, both theoretically (mostly using RRKM, Rice–Ramsperger–Kassel–Marcus) and experimentally, with the reaction Si₂H₆ → SiH₂ + SiH₄ as the most common^53,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,⁶¹ Martin et al.,⁶³ Mick et al.,⁶⁵ Roenigk et al.⁶⁷). 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 H₂Si–SiH₂, H₃Si–SiH, or H₃Si–H₂.^53,54,57 Disilane can even undergo double dehydrogenation to form Si₂H₂, however to our knowledge this has only been reported as the main product when undergoing photolysis at 193 nm.⁶⁹ Yoshida et al.⁵⁷ notes that the transition state for Si₂H₆ → SiH₂ + SiH₄ is 8.48 kcal mol⁻¹ lower than the transition state for Si₂H₆ → H₃SiSiH + H₂, at 43.38 kcal mol⁻¹ compared to 51.86 kcal mol⁻¹.

The spectroscopy of SiH₂ has been used to monitor the SiH₂ + SiH₄ → Si₂H₆ and Si₂H₆ → SiH₂ + SiH₄ reactions (see Fig. 1) and measure the corresponding rate constants and Arrhenius parameters by spectroscopically tracking electronic (ÖX̃) transitions of SiH₂.^39,46,65,70 In these studies, the reconstructions of the amount of SiH₂ 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 SiH₂ molecules after dissociation was made, which could potentially hamper the count of the SiH₂ molecules and affect the reaction rates estimated. A similar experimental technique was used in Hertl and Jolly⁷¹ to monitor SiH₂ in SiH₄ plasma.

It is the second, dissociation, part of the reaction shown in Fig. 1 (Si₂H₆ → SiH₂ + SiH₄) 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 SiH₂ can bear strong non-LTE character, very different from the LTE spectrum of an LTE SiH₂ 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 SiH₂ under conditions similar to dissociation processes expected in these experiments are modelled.

Recently we have computed an accurate ro-vibrational line list for SiH₂, named CATS.⁷² 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.⁷³ 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.

The study by Clark et al.⁷² 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 SiH₂ produced from dissociation of Si₂H₆ through different reaction channels. Using a simplified 1D Harmonic oscillator wavefunctions (see Pastorek et al.⁶) and more sophisticated 3D vibrational CATS wavefunctions from accurate variational calculations, the non-LTE ro-vibrational populations of SiH₂ 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 Si₂H₆ 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)-F12b^76,77 employing the program MOLPRO2015.⁷⁸

Theoretically, the non-LTE properties of dissociating molecules were studied by Band and Freed.⁷⁹ In the present work we use general approach of Berry⁸⁰ and Band and Freed,⁷⁹ which assume no significant structural changes between the reactant and product nuclear configuration, along with the slow vibrational relaxation of the product¹³ to investigate non-LTE spectroscopic signatures of SiH₂ produced from dissociation of Si₂H₆. 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.¹³

We also investigate possible non-LTE impact on the electronic Ã(0,2,0)–X̃(0,0,0) spectrum of SiH₂. This is a favorite spectroscopic system for the detection of SiH₂ due to the large Franck–Condon factor and the availability of suitable laser.^{39,70,81–91}

The non-LTE absorption spectra of SiH₂ are simulated using the (non-LTE) ExoCross program,⁹² 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

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 SiH₂ corresponding to different dissociation routes. In Section 4.2 we calculate the populations and subsequent non-LTE ro-vibrational spectra of SiH₂ using the full 3D wavefunctions and describe the new TROVE methodology. A non-LTE electronic Ã(0,2,0)–X̃(0,0,0) spectrum of SiH₂ is presented in Section 4.4. Conclusions are offered in Section 5.

2 Geometry optimisation and reaction topology of Si₂H₆

2.1 Disilane isomers

In order to better understand the reaction process of breaking Si₂H₆, the topology of Si₂H₆ has been investigated by performing a structural analysis of Si₂H₆ 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 Si₂H₆ helps to indicate how likely local minima or transitions states were to be formed based on corresponding topology. These properties of Si₂H₆ were obtained using the geometry optimization and reaction path finder implemented in MOLPRO2015⁷⁸ using the explicitly correlated coupled cluster method CCSD(T)-F12b^76,77 with the F12-optimized correlation consistent basis set, VTZ-F12⁹⁷ in the frozen core approximation. The calculations employed the diagonal fixed amplitude ansatz 3C(FIX)⁹⁸ and a Slater geminal exponent value of β = 1.0a₀⁻¹.⁹⁹ The auxiliary basis sets were chosen to be the resolution of the identity OptRI¹⁰⁰ basis and the aug-cc-pV5Z/JKFIT¹⁰¹ and cc-awCV5Z/MP2FIT¹⁰² basis sets for density fitting. In the following this level of theory will be referenced to as VTZ/CCSD(T)-F12b.

We shall refer to different disilane isomers as dGM, dLM and dTS to distinguish them from the SiH₂ fragments GM, LM and TS as discussed below.

The global minimum of SiH₂ (dGM) has a symmetrical, staggered D_3d structure. The closest local minimum (dLM) has an inverted, C_3v structure. The lowest transition state, which will be referred to as dTS (TS1 in Becerra et al.⁴⁶ and TS2 in Tonokura et al.¹⁰³), has also inverted structure, just a few kJ mol⁻¹ above dLM.⁴⁶ 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 r_Si–H = 1.492 Å, r_Si–Si = 2.331 Å and bond angles β_∠HSiSi = 110.3° and α_∠HSiSi = 108.6°.¹⁰⁴ The structure of the deuterated isotopologue Si₂H₅D was reported as r_Si–H = 1.4874(17) Å, r_Si–Si = 2.3317(15) Å, β_∠HSiSi = 110.66(16)°.¹⁰⁵

Table 1 An overview of the Si₂H₆ molecules discussed in the paper. The subscript L denotes atoms to the LHS of the Si–Si bond in the figures shown, and the subscript R denotes atoms to the RHS

Isomer	LM	TS	GM	Becerra et al.^46
Molecule				LM1	TS1	GM
Si–Si (Å)	2.467	2.489	2.341	2.470	2.485	2.346
SiR–HR (Å)	1.532	1.499	1.482	1.525	1.493	1.479
SiL–HL (Å)	1.480	1.480	1.482	1.478	1.478	1.479
∠HRSiSi (°)	55.4	64.9	110.2	—	—	110.24
∠HLSiSi (°)	110.7	110.0	110.2	110.7	110.2	110.24
∠HRSiHR (°)	90.9	100.6	108.6	91.4	102.6	—
∠HLSiHL (°)	108.3	109.7	108.6	108.2	108.8	—
ω b (cm^−1)	940	940	940
ω s(L) (cm^−1)	2241.69	2245.36	2240.10
ω s(R) (cm^−1)	2039.12	2170.66	2240.10

---

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⁻¹ and 4.1 kJ mol⁻¹ 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.⁴⁶ and Sakai and Nakamura,⁶⁰ albeit both vary for the dGM structure by 36 kJ mol⁻¹.

Fig. 2 Potential energy curve for the reaction path from disilane global minimum (dGM), left, through the transition state (dTS) onto the local minimum (dLM), right. A zoom of dLM to dTS path is inset. The dLM is 344.22 cm⁻¹ lower than the energy of the dTS.

Table 2 The total and relative energies for the three disilane structure. Literature values from Becerra et al.^c 46 and Sakai and Nakamura⁶⁰

	Energy, Eh	Rel. energy,^d kJ mol^−1
a TS2 in Sakai and Nakamura^60 is the same as TS1 in Becerra et al.^46 b Complex 2 in Sakai and Nakamura^60 is the same as LM1 in Becerra et al.^46 c Relative to SiH2 + SiH4, kJ mol^−1. d Relative to dTS.
	dGM
Calculated	−581.728134	−178.00
Becerra	−236.0	−214.50
Saki	−581.53962	−214.68715
	dTS
Calculated	−581.660338	0.00
Becerra (TS1)^a	−48.5	0.00
Saki (TS2)^a	−581.45785	0.00
	dLM
Calculated	−581.661906	−4.12
Becerra (LM1)^b	−51.3	−3.00
Saki (compl. 2)^b	−581.45837	−1.3652601

---

The dLM isomer of Si₂H₆ has a shallow potential with a very low barrier to dTS of 344 cm⁻¹. It can be also recovered using lower levels of theory, for example, using MP2/6-311G(d,p)⁴⁶ and even with the UFF force fields implemented in Avogadro 1.2.0¹⁰⁶via the steepest descent method and 4 steps per update.

2.2 The silylene fragments

With the aim to give more quantitative information on the structural and dynamical properties of five SiH₂ fragments from Si₂H₆, the fragments are described as follows. GM is an SiH₂ fragment from the global minimum structure (dGM); LM-L is an SiH₂ fragment from the left hand side (LHS) of the local minimum structure (dLM) with the Si–H₃ umbrella group pointing outside; LM-R is an SiH₂ fragment from the right hand side (RHS) of the dLM structure with Si–H₃ pointing inside; TS-L is an SiH₂ fragment from the LHS of the transition state structure (dTS), Si–H₃ umbrella group points outside; and TS-R is an SiH₂ fragment from the RHS of the dTS structure, Si–H₃ 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 D_3d whereas the symmetry for the dLM and dTS are both C_3v. 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⁻¹.

Table 3 The calculated and available literature vibrational frequencies in cm⁻¹ for the three Si₂H₆ molecules

dGM, D3d				dTS, C3v				dLM, C3v
Calculated		Literature^46		Calculated		Literature^46,67		Calculated
Mode	Freq.	Mode	Freq.	Mode	Freq.	Mode	Freq.	Mode	Freq.
a Confirmed by the experimental gas phase Raman spectrum of Durig and Church.^59 b 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).
E u	2247.78	E u	2179	E	2252.06	A	2181	E	2244.30
A 1g	2241.03	A 1g	2163^a	A	2238.66	E	2169	A	2239.08
E g	2239.28	E u	2155	E	2176.47	A	2105	A	2050.89
A 1u	2232.30	A 1u	2154	A	2161.85	A	2087	E	2027.35
E u	962.35	E g	941^a	E	995.37	A	1585	E	1089.15
E g	948.37	E u	940	E	966.95	A	960	E	978.33
A 1g	931.91	A 1g	920^a	A	902.23	A	949	A	956.33
A 1u	859.02	A 1u	844	E	570.16	A	927	A	667.30
E g	636.87	E g	628^a	A	353.88	A & E^b	925	E	618.92
A 1g	438.76	A 1g	432^a	E	319.45	E	145	A	432.13
E u	372.85	E u	379	A	316.43	A	128	A	389.31
A 1u	136.56	A 1u	128	A	−551.59	E	93	E	277.91

---

3 Modelling the non-LTE populations of SiH₂

We now consider the decomposition reaction Si₂H₆ → SiH₂ + SiH₄ shown in Fig. 1 and model the vibrational populations of the product SiH₂ 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 SiH₂ using these populations to demonstrate their impact on the spectral shape of dissociated species.

In our description of the non-LTE vibrational population of the dissociated molecule we follow the polyatomic Franck–Condon type approximation by Band and Freed⁷⁹ and Berry⁸⁰ 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 SiH₂) to an initial structure of Si₂H₆ through the dissociation process, we assume that the dissociation happens instantaneously, i.e. the initial configuration of the product SiH₂ corresponds to the structural parameters (bond lengths Si–H and inter-bond angles ∠HSiH) of SiH₂ as a fragment of Si₂H₆, 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 SiH₂ is assumed to be r_Si–H = 1.482 Å, α_∠HSiH = 108.6°. Naturally, this is a very deformed geometry comparing to the equilibrium structure of the gas phase SiH₂, r_e = 1.5137 ± 0.003 Å and α_e = 92.04 ± 0.05°.¹⁰⁷ After being dissociated, in relaxing to be a free molecule, the fragment SiH₂ 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.

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.¹³). The vibrational population however is not in the thermal equilibrium and therefore no sensible vibrational temperature could be associated with the corresponding population.

The non-LTE population N_J,k,ν(T) of a ro-vibrational state |J, k, ν〉 is then given by:⁶

where N^vib_ν is a non-LTE vibrational population, T_rot 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 (ν₁,ν₂,ν₃) to describe vibrational states of a triatomic molecule; g^ns_i is the nuclear spin degeneracy; T is the temperature, c₂ is the second radiation constant. In eqn (1), Ẽ^rot_J,k is the rotational part of the ro-vibrational energy term value approximated as

Ẽ_J,k,ν = Ẽ^vib_ν + Ẽ^rot_J,k, (2)

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 Q_nLTE in eqn (1) is given byOur aim is to calculate the vibrational populations N_v ≡ N^vib_ν of SiH₂ as produced by instantaneous (vertical) dissociation from three structures, dGM, dLM and dTS. In case of dLM and dTS, the SiH₂ 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.

Table 4 Five SiH₂ fragments considered for the non-LTE analysis

SiH2 fragment name	Formation & structure
GM	From dGM, two equivalent SiH3
TS-R	From dTS, SiH3 pointing inside
TS-L	From dTS, SiH3 pointing outside
LM-R	From dLM, SiH3 pointing inside
LM-L	From dLM, SiH3 pointing outside

---

Let us assume that the Si₂H₆ 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 ≡ (ν₁,ν₃) with some transition probability giving rise to the vibrational population N^vib_ν. On top of that we also assume a full separation of the stretching Si–H and bending ∠HSiH modes inside Si₂H₆ in its ground vibrational state. Possible consequences of deviation from these approximations are discussed below.

Under the assumptions made we define the vibrational population of the SiH₂ fragment as a Franck–Condon factor for a vertical transition from the ground vibrational state of disilane Si₂H₆ to SiH₂ + SiH₄ with the gas phase (g.ph.) SiH₂ transferred to some vibrational state |(g.ph.)〉 = |ν₁,ν₂,ν₃〉. In the approximation of the full separation of the fragment SiH₂ from the rest of Si₂H₆, the population of SiH₂ can be represented as an overlap between the ground state wavefunction |ν′′ = 0(fragment)〉 of a fragment SiH₂ and that of the corresponding vibrational state |ν′ = (g.ph.)〉 of the gas phase SiH₂ as given by

N = |〈0(fragment)|ν(g.ph.)〉|². (4)

Here the spacial wavefunction of the SiH₂ fragment (i.e. a combination of two adjacent Si–H bonds in disilane with an angle between them forming the dissociating SiH₂) is projected on vibrational eigenfunctions of the gas phase SiH₂ to give the corresponding populations of SiH₂. The highest populated energy level will have the largest overlap between these wavefunctions.

The calculated temperature dependent populations N_J,k,ν(T) in eqn (1) can be then combined with a molecular line list for SiH₂ to simulate absorption or emission spectra of this molecule under the non-LTE conditions as defined by N_J,k,ν(T) in eqn (1). Here we use the ExoMol line list CATS by Clark et al.⁷² as provided by ExoMol (http://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 elsewhere⁷³). A non-LTE spectrum of SiH₂ for given T and P is then calculated using CATS' Einstein A coefficients with the ExoCross program,⁹² 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).

4 Computing vibrational populations of SiH₂

Two approaches were used for the calculation of the population densities of the fragment SiH₂. One approach – named the decoupled 1D approach – is where the 3D wavefunctions of the fragment as well as of free SiH₂ 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.⁷⁴ Both approaches are presented in the following in order to assess and compare the accuracy achieved.

4.1 The 1D approach

A vibrational state |ν〉 = |ν₁,ν₂,ν₃〉 of SiH₂ is characterized by the three (normal mode) quantum numbers ν₁, ν₂ and ν₃ corresponding to the two stretching modes (ν₁ and ν₃) and one bending mode (ν₂) of the SiH₂ molecule. The 1D approach considers the stretching Si–H₁, Si–H₂ and bending ∠HSiH modes, both of the molecular fragment and gas phase SiH₂ molecules, as fully independent and described by one dimensional (1D) wavefunctions under the harmonic approximation, as given by:Here x is a dimensionless coordinate describing either the stretching r = r_Si–H or bending α = α_{∠H–Si–H} coordinate as follows:withand ω_s = ω_Si–H, ω_b = ω_∠HSiH. In eqn (5)H(x) is a Hermite polynomial and C_v is the corresponding normalization constant.

The constants a_s and a_b correspond to inverse masses of the vibrational part of a free three-atomic molecule expressed in terms of the internal coordinates r₁, r₂, α (see, e.g. Sutcliffe and Tennyson,¹⁰⁸ Yurchenko et al.¹⁰⁹).

A 1D population for the corresponding mode of the gas phase SiH₂ molecule is given by eqn (4) with |ν〉 = Ψ_v(x). The different disilane fragments have different structural parameters r_e, α_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 SiH₂ according with eqn (4).

For each of the three modes (two stretching and one bending), 1D wavefunctions of the gas phase SiH₂ for 30 vibrational states from ν₀ up to ν₂₉ were calculated. These 1D wavefunctions were then numerically integrated with the corresponding ground state 1D wavefunctions |ν = 0〉 of the fragment in question.

The total vibrational population N_ν1,ν2,ν3 in this approximation is then given by a product

where and are obtained using the stretching harmonic oscillators wavefunctions |ν₁〉 and |ν₃〉, while is obtained using the corresponding bending harmonic oscillator wavefunction |ν₂〉. The independent treatment of the two stretching populations is partly justified by the local mode character of the vibrational degrees of freedom of SiH₂ due to the near 90° bond angle (see, e.g. Jensen¹¹⁰ and Clark et al.⁷²). The asymmetric vibrational modes of SiH₂ (B₂ in C_2v) 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 (A₁), only symmetric states of SiH₂ 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 A₁) 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 B₂-type (0,0,1) vibratitonal state is assumed to be zero. The populations N_ν1,ν2,ν3 are pre-calculated for each vibrational state (ν₁,ν₂,ν₃) of SiH₂ 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.†

4.1.1 1D populations and spectra. 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 SiH₂.⁷² 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 SiH₂. 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.

Fig. 3 The 1D bending mode harmonic wavefunctions for the five SiH₂ fragments described in Table 1 as a function of bond angle in degrees. Ground state fragment (black) compared with the LTE |ν = 0〉 (blue), |ν = 1〉 (green), |ν = 2〉 (red).

Fig. 4 The 1D stretching mode harmonic wavefunctions for the five SiH₂ fragments described in Table 1 as a function of bond length in Angstrom. Ground state fragment (black) compared with the LTE |ν = 0〉 (blue), |ν = 1〉 (green), |ν = 2〉 (red).

Fig. 5 Stretching and bending mode 1D populations N_v for the five SiH₂ fragments from ν = 0–ν = 10.

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.

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.

Fig. 6 and 7 show the 1000 cm⁻¹ and 2000 cm⁻¹ (10 μm and 5 μm) bands for the SiH₂ 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⁻¹ 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.

Fig. 6 Non-LTE spectra of SiH₂ at T_rot = 296 K in the 1000 cm⁻¹ (10 μm) region corresponding to five vibrational populations GM, TS-L, TS-R, LM-L and LM-R (upper displays of each figure) and compared to the same spectra simulated using the LTE population at T = 296 K. The non-LTE populations were obtained using the 1D approach (see text). Only the strongest bands are shown. The CATS line list was used.

Fig. 7 Non-LTE spectra of SiH₂ at T_rot = 296 K in the 2000 cm⁻¹ (5 μm) region corresponding to five vibrational populations GM, TS-L, TS-R, LM-L and LM-R (upper displays of each figure) and compared to the same spectra simulated using the LTE population at T = 296 K. The non-LTE populations were obtained using the 1D approach (see text). Only the strongest bands are shown. The CATS line list was used.

The plots in Fig. 7 show the 2000 cm⁻¹ 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).

With the equilibrium structures of the TS-R and LM-R fragments being similar to the equilibrium structure of SiH₂, 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⁻¹ and 2000 cm⁻¹ bands.

4.2 The 3D approach for vibrational populations using an accurate variational method

In a full 3D approach, the ground state wavefunction |0,0,0(fragment)〉 represents an SiH₂ fragment of an 18D ground state vibrational wavefunction of Si₂H₆:

ϕ_Si2H6 = |0,0,0(SiH₂)〉 |0,0,0,0,0,0,0,0,0,0,0,0,0,0,0(SiH₄)〉.

Here we assume the approximation that the corresponding three modes (Si–H₁, Si–H₂ and ∠H₁SiH₂) are independent from the rest of the molecule so that all other modes, not relevant for the gas phase SiH₂, can be eliminated (integrated out), including the reaction coordinate and vibrational modes of SiH₄. This is in line with the assumptions used previously by Band and Freed⁷⁹ and Berry.⁸⁰ Apart from this approximation we will treat the SiH₂ 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).

The vibrational populations are then modelled using the Franck–Condon integrals as follows:

N_ν1,ν2,ν3 = |〈0,0,0(fragment)|ν₁, ν₂, ν₃(g.ph.)〉|², (11)

with the sum of all populations over all states totalling 1. In this equation, |ν₁,ν₂,ν₃(g.ph.)〉 is an accurate vibrational wavefunction of a gas phase SiH₂ molecule, obtained by solving the vibrational Schrödinger equation with an accurate PES. We use the TROVE variational program and the refined PES of SiH₂ by Clark et al.⁷² to generate |ν₁,ν₂,ν₃(g.ph.)〉 for all vibrational excitations required. For the details on the TROVE calculations see below and also Clark et al.⁷²

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 SiH₂. 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 products

of the corresponding eigen-coefficients C_i1,i2,i3(fragment) and C_i1,i2,i3(g.ph.), obtained variationally in independent calculations using a new implementation in TROVE.

TROVE uses optimized non-standard vibrational basis sets, generated numerically by solving 1D Schrödinger equations for realistic 1D potentials.⁷⁴ 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 SiH₂ 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.

4.2.1 PESs of SiH₂ fragments. For our new 3D populations corresponding to dissociations from the fragments, five PESs were generated as follows. We assume that PESs of a dissociating Si₂H₆ molecule can be approximated as a sum of two independent fragments:

V_Si2H6 = V_SiH2(r₁, r₂, α) + V_rest,

where the individual stretching and bending modes of SiH₂ fragments are fully separable:

V_SiH2(r₁, r₂, α) = f_str(r₁) + f_str(r₂) + f_bnd(α). (13)

The stretching part of the potential is given by a Morse-like expansionwhile the bending part is a Taylor-type expansion in terms of the displacement from the corresponding equilibrium value:The expansion constants f_i and g_i representing fragments' potential energies V_SiH2(r₁,r₂,α) 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 Si₂H₆ (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.

Table 5 The ab initio potential parameters defining the 1D potential of different structures used in TROVE calculations. The units are cm⁻¹, Å and radians, unless specified

Parameter	GM	TS-L	TS-R	LM-L	LM- R
r e, Å	1.48218	1.47971	1.49888	1.47982	1.51611
α e, deg	108.701	108.938	103.268	108.264	92.287
a, Å^−1	1.29065	1.29110	1.22364	1.27911	1.27050
g 1	0	0	0	0	−8132.4642
g 2	80964.377	90039.635	101613.84	93608.413	168296.56
g 3	−117331.83	−126539.98	−159054.10	−129420.43	−299843.15
g 4	139212.98	168090.27	182459.45	305247.60	324286.87
f 1	0	0	0	0	−793.07685
f 2	43617.985	43796.459	45462.423	22217.055	19382.230
f 3	−6594.4391	−6937.7586	−8917.729	−3824.8354	−2298.2220
f 4	2809.6655	2714.6696	3750.0966	1407.7267	2233.5205

---

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 SiH₂ 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.

Fig. 8 1D potential energy cuts representing different stretching Si–H (x-axis = bond length) and bending (x-axis = bond angle) H–Si–H modes of different fragments of three Si₂H₆ isomers, dGM, dLM and dTS, compared to the corresponding cuts of the gas phase SiH₂ species.

4.2.2 Vibrational calculations. The vibrational wavefunctions of SiH₂ were computed using the variational nuclear motion program TROVE with the same setup as in Clark et al.⁷² 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 SiH₂ equilibrium geometry in terms of linearized coordinates.¹¹⁴ The primitive basis set is constructed from 1D mode numerical basis functions using the Numerov–Cooley approach^115,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 C_2v(M) symmetry group¹¹⁴ using an optimized symmetrization procedure¹¹² 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-off

P= 2(ν₁ + ν₃) + ν₂ ≤ 24, (16)

where ν₁ and ν₃ are the stretching and ν₂ is the bending quantum numbers, with the maximal excitations 12, 12 and 24 respectively.

For the gas phase SiH₂ calculations we employed the empirically refined PES by Clark et al.⁷² For the five SiH₂ 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).

4.2.3 Vibrational populations of SiH₂. The 3D vibrational populations N_ν1,ν2,ν3 in eqn (11) were computed using the corresponding TROVE eigenfunctions viaeqn (12).

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 SiH₂ 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 SiH₂ 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.

Fig. 9 3D Bending primitive wavefunctions of SiH₂ compared to the ground state 3D bending wavefunctions for Fragment TS-L. Calculated with TROVE.

Fig. 10 The 3D non-LTE vibrational populations of SiH₂ for the dissociation from TS-L (blue circles) as a function of the vibrational term values compared to Boltzmann distribution of SiH₂ vibrational states at T = 1000 K.

4.3 Non-LTE intensity simulations

A non-LTE absorption line intensity I_fi (cm per molecule) can be calculated aswhere A_fi is the Einstein A coefficient (s⁻¹), [small nu, Greek, tilde]_fi is the transition wavenumber (cm⁻¹), Q_nLTE(T) is the non-LTE partition function defined in eqn (3).

In order to simulate absorption spectra of the gas phase SiH₂ assuming a non-LTE vibrational populations N_ν1,ν2,ν3, we use the line positions and Einstein A coefficients from the ExoMol CATS line list⁷² employing the ExoCross program.⁹² The ExoMol line lists are formatted as two files, a States file and a Transition file. It is described extensively elsewhere⁷³ and in this paper we shall only discuss how the ExoMol format has been adapted for use in non-LTE situations.

To adapt the CATS States file⁷² 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 SiH₂. The ‘density’ column is specific for the calculation of line lists for non-LTE molecules is not routinely included into the ExoMol States files.

Table 6 Extract from the modified.states file of the GM 1D non-LTE line list

i: state counting number; Ẽ_i: state energy in cm⁻¹; g_i: state degeneracy; J: total angular momentum quantum number; Γ_tot: overall symmetry of state in C_2v (M); ν₁–ν₃: vibrational (normal mode) quantum numbers; Γ_vib: vibrational symmetry in C_2v (M); K_a: asymmetric top quantum number; K_c: asymmetric top quantum number; Γ_rot: rotational symmetry in C_2v(M); C_i: largest coefficient used in the TROVE assignment; n₁–n₃: vibrational (TROVE) quantum numbers; i_vib: vibrational state counting number; N^vib: population density for each i based on ν₁–ν₃. This column is not produced by the ExoMol format. It must be calculated separately.

---

This column is read by ExoCross and used to give the population weighting to each line intensity. Evaluation of the non-LTE population N_J,k,ν(T) given for each ro-vibrational state as in eqn (1) is based on the knowledge of the corresponding vibrational state (ν₁,ν₂,ν₃) as well as the rotational energy contribution Ẽ^rot_J,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 Ẽ_v^vib from the total energy according with eqn (2). All vibrational quantum numbers (ν₁,ν₂,ν₃) are not required, only a vibrational index indicating the vibrational state in question. In our model, all asymmetric vibrational states (B₂) are not populated due to the zero overlap with the ground state of the A₁ symmetry in eqn (11), as part of the completely vertical Franck–Condon approximation.

4.3.1 Using the non-LTE populations from the 3D approach. Fig. 11 shows non-LTE spectra of SiH₂ in the two main spectroscopic regions, 1000 and 2000 cm⁻¹ 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 T_rot = 296 K. These non-LTE spectra are compared to the corresponding LTE spectra of SiH₂ at T = 296 K, which comprise mainly two fundamental bands, (0,1,0)←(0,0,0) (1000 cm⁻¹ region) and (0,0,1)←(0,0,0) (2000 cm⁻¹ 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⁻¹ band), (1,1,0)←(0,1,0) and (1,1,0)←(0,2,0) (2000 cm⁻¹ 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⁻¹ region, with the difference of about 15 cm⁻¹. 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 r_SiH and α_{∠H–Si–H} of the LM-R structure are very similar to the corresponding equilibrium parameters of the gas phase SiH₂.

Fig. 11 Non-LTE spectra of SiH₂ at T_rot = 296 K corresponding to five vibrational populations GM, TS-L, TS-R, LM-L and LM-R (upper displays of each figure) and compared to the same spectra simulated using the LTE population at T = 296 K. The non-LTE populations were obtained using the 3D TROVE approach. Only three or four strongest bands are shown. The CATS line list was used.

The stark differences of the non-LTE spectra offer an ability for experiment to distinguish between SiH₂ molecules produced from fragmenting Si₂H₆, and even to indicate dissociation channels involved.

4.4 ÖX̃ spectrum

The visible electronic band ÖX̃ of SiH₂ has often been used to study different reactions involving leading to silylene.⁹¹ Here we used the program RENNER¹¹⁷ to simulate a non-LTE electronic spectrum of SiH₂ with the spectroscopic model by Yurchenko et al.¹¹⁸ for the ùB₁–X̃¹A₁ 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).

In the RENNER calculations, the size of the basis set originally used in Yurchenko et al.¹¹⁸ 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̃¹A₁ electronic state basis set included N_A = 12 stretching functions of the A₁ symmetry and N_B = 10 stretching functions of the B₂ symmetry. For the ùB₁ state N_A = 10 and N_B = 8 stretching functions were used. These stretching functions were constructed from the Morse oscillator functions |n₁〉|n₃〉 with n₁ + n₃ ≤ N_stretch = 12.

A rovibronic line list for the ùB₁–X̃¹A₁ of SiH₂ was generated covering the rotational excitations up to J_max = 15 with the lower state energies (X̃) truncated at hc·25 000 cm⁻¹ and the upper state energies truncated at hc·28 000 cm⁻¹.

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 SiH₂ in the region of the band Ã(0,2,0)←X̃(0,0,0), assuming the rotational temperature T_rot = 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 1₀₁←1₁₀ belonging to this band used in a number of experimental studies involving SiH₂ as a reaction product^71,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 SiH₂ is required to estimate the number density of SiH₂ 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 SiH₂ as a reaction product can vary significantly under non-LTE depending on the reaction pathway and impact the experimental rates.

Fig. 12 An ùB₁–X̃¹A₁ spectra of SiH₂, LTE (T = 500 K) and non-LTE T_rot = 500 K using the 1D model for Fragment TS-L. The energies and Einstein coefficients are generated using the RENNER program with the spectroscopic model from Yurchenko et al.¹¹⁸ A Gaussian line profile of HWHM = 5 cm⁻¹. The popular ro-vibronic line (0,0,0)1₀₁←(0,0,0)1₁₀ used in measurements of reaction rates of SiH_n species^71,89,90 is shown.

5 Conclusion

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 SiH₂ calculated with 3D wavefunctions and TROVE, and compared them to non-LTE spectra of SiH₂ modelled using a 1D harmonic approach, and the LTE line list calculated previously by ExoMol.

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 SiH₂ corresponding to dissociation of disilane from different sides of the three disilane isomer were computed. We have shown that the non-LTE spectra of SiH₂ are different in most cases. This is important as the spectrum of SiH₂ is used to monitor the quantity present in a reaction as a means to track the progress of SiH₂ + SiH₄ → Si₂H₆ and Si₂H₆ → SiH₂ + SiH₄ when calculating the corresponding rate constant. If Si₂H₆ is decomposing at a rate slower than it is being formed, then tracking the quantity of SiH₂ 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 SiH₂ and Si₂H₆.

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 SiH₂ fragment during the instantaneous dissociation and is fully decoupled from the rest of the Si₂H₆ molecule, i.e. can be described by a 3D wavefunction in its lowest, relaxed vibrational configuration and has the same structural parameters as the Si₂H₆ molecule.

We have shown that the non-LTE spectra of SiH₂ 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 SiH₄, which has not been explored here.

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,⁶ recently, which are planning to explore in the future.

The methods described here can be used model the intensity distribution of the reaction products and to ascertain from what molecule the SiH₂ dissociated from. The equilibrium structured parameters (bond lengths and angles) can be treated as effective parameters to be adjusted to reproduce the experimental spectra.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

This work was supported by UK research councils EPSRC, under grant EP/N509577/1 with COVID extension, and STFC, under grant ST/R000476/1. This work made extensive use of the STFC DiRAC HPC facility supported by BIS National E-infrastructure capital grant ST/J005673/1 and STFC grants ST/H008586/1 and ST/K00333X/1. We thank the European Research Council (ERC) under the European Unions Horizon 2020 research and innovation programme through Advance Grant number 883830. We also thank Thomas Mellor for help with the variational model.

References

1. H. W. Hermann and S. R. Leone, J. Chem. Phys., 1982, 76, 4766–4774.
2. M. Ferus, P. Kubelík, K. Kawaguchi, K. Dryahina, P. Španěl and S. Civiš, J. Phys. Chem. A, 2011, 115, 1885–1899.
3. M. Ferus, R. Michalčíková, V. Shestivská, J. Šponer, J. E. Šponer and S. Civiš, J. Phys. Chem. A, 2014, 118, 719–736.
4. M. Ferus, P. Kubelík, A. Knížek, A. Pastorek, J. Sutherland and S. Civiš, Sci. Rep., 2017, 7, 6275.
5. K. Prozument, J. H. Baraban, P. B. Changala, G. B. Park, R. G. Shaver, J. S. Muenter, S. J. Klippenstein, V. Y. Chernyak and R. W. Field, Proc. Natl. Acad. Sci. U. S. A., 2020, 117, 146–151.
6. A. Pastorek, S. Civiš, V. H. J. Clark, S. N. Yurchenko and M. Ferus, J. Quant. Spectrosc. Radiat. Transf., 2021, 262, 107521.
7. P. G. Wenthold, D. A. Hrovat, W. T. Borden and W. C. Lineberger, Science, 1996, 272, 1456–1459.
8. W. H. Green Jr, C. B. Moore and W. F. Polik, Annu. Rev. Phys. Chem., 1992, 43, 591–626.
9. E. A. Michael, C. J. Keoshian, S. K. Anderson and R. J. Saykally, J. Mol. Spectrosc., 2001, 208, 219–223.
10. E. Michael, C. Keoshian, D. Wagner, S. Anderson and R. Saykally, Chem. Phys. Lett., 2001, 338, 277–284.
11. G. Bazalgette Courrèges-Lacoste, J. Sprengers, J. Bulthuis, S. Stolte, T. Motylewski and H. Linnartz, Chem. Phys. Lett., 2001, 335, 209–214.
12. D. Zhao, K. D. Doney and H. Linnartz, J. Mol. Spectrosc., 2014, 296, 1–8.
13. E. Dudás, N. Suas-David, S. Brahmachary, V. Kulkarni, A. Benidar, S. Kassi, C. Charles and R. Georges, J. Chem. Phys., 2020, 152, 134201.
14. E. A. Milne, Handbuch der Astrophysik, Springer, 1930, vol. 3, ch. 2, p. 65.
15. A. Curtis and R. Goody, Proc. Math. Phys. Eng., 1956, 236, 193–206.
16. J. Houghton, Q. J. R. Meteorol. Soc., 1969, 95, 1–20.
17. W. R. Kuhn and J. London, J. Atmos. Sci., 1969, 26, 189–204.
18. R. E. Dickinson, J. Atmos. Sci., 1972, 29, 1531–1556.
19. J. B. Kumer and T. C. James, J. Geophys. Res., 1974, 79, 638–648.
20. G. Shved, Soviet Astron., 1975, 18, 499–504.
21. M. López-Puertas, R. Rodrigo, J. López-Moreno and F. Taylor, J. Atmos. Terr. Phys., 1986, 48, 749–764.
22. M. López-Puertas, R. Rodrigo, A. Molina and F. Taylor, J. Atmos. Terr. Phys., 1986, 48, 729–748.
23. P. P. Wintersteiner, R. H. Picard, R. D. Sharma, J. R. Winick and R. A. Joseph, J. Geophys. Res. - Atmos., 1992, 97, 18083–18117.
24. A. P. Williams, PhD thesis, University of Oxford, 1971.
25. F. van der Tak, J. Black, F. Schöier, D. Jansen and E. van Dishoeck, Astron. Astrophys., 2007, 468, 627–635.
26. B. Funke, M. López-Puertas, M. García-Comas, M. Kaufmann, M. Höpfner and G. Stiller, J. Quant. Spectrosc. Radiat. Transf., 2012, 113, 1771–1817.
27. E. Pannier and C. O. Laux, J. Quant. Spectrosc. Radiat. Transf., 2019, 222–223, 12–25.
28. M. López-Puertas and F. W. Taylor, Non-LTE radiative transfer in the Atmosphere, World Scientific, 2001, vol. 3.
29. M. Asplund, Annu. Rev. Astron. Astrophys., 2005, 43, 481–530.
30. S. Chandrasekhar, Radiative Transfer, Dover Publications, New York, 1960.
31. P. Feautrier, Comptes Rendus Acad. Sci. Paris, 1964, 258, 3189.
32. G. Rybicki, J. Quant. Spectrosc. Radiat. Transf., 1971, 11, 589–595.
33. R. G. Athay, Radiation transport in spectral lines, D. Reidel Publishing Co. (Geophysics and Astrophysics Monographs), Dordrecht, Holland, 1972, vol. 1.
34. D. Mihalas, Stellar Atmospheres, WH Freeman and Co, San Francisco, 1978.
35. E. Simonneau and L. Crivellari, Astrophys. J., 1993, 409, 830.
36. H. Chen, Z. Zhan, A. Youngblood, E. T. Wolf, A. D. Feinstein and D. E. Horton, Nat. Astron., 2021, 5, 298–310.
37. P. F. Goldsmith and W. D. Langer, Astrophys. J., 1999, 517, 209–225.
38. F. Lique, F. van der Tak, J. Klos, J. Bulthuis and M. Alexander, Astron. Astrophys., 2009, 493, 557–563.
39. J. Jasinski and J. O. Chu, J. Chem. Phys., 1988, 88, 1678–1687.
40. J. R. Doyle, D. A. Doughty and A. Gallagher, J. Appl. Phys., 1990, 68, 4375–4384.
41. J. Jasinski, B. S. Meyerson and B. A. Scott, Annu. Rev. Phys. Chem., 1987, 38, 109–140.
42. M. Ernst, A. Sax and J. Kalcher, Chem. Phys. Lett., 1993, 216, 189–193.
43. D. M. Goldhaber and A. L. Betz, Astrophys. J., 1984, 279, L55–L58.
44. R. Kaiser and Y. Osamura, Astron. Astrophys., 2005, 432, 559–566.
45. P. M. Agrawal, D. L. Thompson and L. M. Raff, J. Chem. Phys., 1988, 89, 741–750.
46. R. Becerra, H. M. Frey, B. P. Mason, R. Walsh and M. S. Gordon, J. Chem. Soc., Faraday Trans., 1995, 91, 2723–2732.
47. C. Puzzarini and P. R. Taylor, J. Chem. Phys., 2005, 122, 054315.
48. J. Church, J. Durig, T. A. Mohamed and A. B. Mohamad, Spectrochim. Acta, Part A, 1994, 50, 639–660.
49. H. Sun, Macromolecules, 1995, 28, 701–712.
50. A. Márquez, J. F. Sanz, M. Gelizé and A. Dargelos, Chem. Phys., 1991, 149, 311–318.
51. A. Romero, M. Kiwi and R. Ramrez, Phys. Status Solidi B, 2002, 230, 391–395.
52. S. G. Cho, O. K. Rim and G. Park, J. Comput. Chem., 1997, 18, 1523–1533.
53. P. M. Agrawal, D. L. Thompson and L. M. Raff, J. Chem. Phys., 1990, 92, 1069–1082.
54. S.-W. Hu, Y. Wang, X.-Y. Wang, T.-W. Chu and X.-Q. Liu, J. Phys. Chem. A, 2003, 2954–2963.
55. M. M. Law, J. T. Fraser-Smith and C. U. Perotto, Phys. Chem. Chem. Phys., 2012, 14, 6922–6936.
56. W. A. Eger, A. Genest and N. Rösch, Chem. – Eur. J., 2012, 18, 9106–9116.
57. K. Yoshida, K. Matsumoto, T. Oguchi, K. Tonokura and M. Koshi, J. Phys. Chem., 2006, 110, 4726–4731.
58. L. Andrews and X. Wang, J. Phys. Chem. A, 2002, 106, 7696–7702.
59. J. R. Durig and J. S. Church, J. Chem. Phys., 1980, 73, 4784–4797.
60. S. Sakai and M. Nakamura, J. Phys. Chem., 1993, 97, 4960–4965.
61. M. Bowrey and J. Purnell, Proc. Math. Phys. Eng. Sci., 1971, 321, 341–359.
62. G. Olbrich, P. Potzinger, B. Reimann and R. Walsh, Organometallics, 1984, 3, 1267–1272.
63. J. G. Martin, H. E. O'Neal and M. A. Ring, Int. J. Chem. Kinet., 1990, 22, 613–631.
64. H. K. Moffat, K. F. Jensen and R. W. Carr, J. Phys. Chem., 1992, 96, 7683–7695.
65. H.-J. Mick, M. W. Markus, P. Roth and V. N. Smirnov, Ber. Bunsenges. Phys. Chem., 1995, 99, 880–890.
66. A. A. Onischuk and V. N. Panfilov, Russ. Chem. Rev., 2001, 70, 321–332.
67. K. F. Roenigk, K. F. Jensen and R. W. Carr, J. Phys. Chem., 1987, 91, 5732–5739.
68. H. Emeléus and C. Reid, J. Chem. Soc., 1939, 1021–1030.
69. N. Tada, K. Tonokura, K. Matsumoto, M. Koshi, A. Miyoshi and H. Matsui, J. Phys. Chem. A, 1999, 103, 322–329.
70. M. W. Markus and P. Roth, J. Quant. Spectrosc. Radiat. Transf., 1994, 52, 783–789.
71. M. Hertl and J. Jolly, J. Phys. D: Appl. Phys., 2000, 33, 381–388.
72. V. H. J. Clark, A. Owens, J. Tennyson and S. N. Yurchenko, J. Quant. Spectrosc. Radiat. Transf., 2020, 246, 106929.
73. J. Tennyson, S. N. Yurchenko, A. F. Al-Refaie, V. H. J. Clark, K. L. Chubb, E. K. Conway, A. Dewan, M. N. Gorman, C. Hill, A. E. Lynas-Gray, T. Mellor, L. K. McKemmish, A. Owens, O. L. Polyansky, M. Semenov, W. Somogyi, G. Tinetti, A. Upadhyay, I. Waldmann, Y. Wang, S. Wright and O. P. Yurchenko, J. Quant. Spectrosc. Radiat. Transf., 2020, 255, 107228.
74. S. N. Yurchenko, W. Thiel and P. Jensen, J. Mol. Spectrosc., 2007, 245, 126–140.
75. S. N. Yurchenko, A. Yachmenev and R. I. Ovsyannikov, J. Chem. Theory Comput., 2017, 13, 4368–4381.
76. T. B. Adler, G. Knizia and H. J. Werner, J. Chem. Phys., 2007, 127, 221106.
77. G. Knizia, T. B. Adler and H.-J. Werner, J. Chem. Phys., 2009, 130, 054104.
78. H.-J. Werner, P. J. Knowles, G. Knizia, F. R. Manby and M. Schütz, Wiley Interdiscip. Rev.: Comput. Mol. Sci., 2012, 2, 242–253.
79. Y. B. Band and K. F. Freed, J. Chem. Phys., 1975, 63, 3382–3397.
80. M. J. Berry, Chem. Phys. Lett., 1974, 27, 73–77.
81. G. Inoue and M. Suzuki, Chem. Phys. Lett., 1984, 105, 641–644.
82. J. M. Jasinski, E. A. Whittaker, G. C. Bjorklund, R. W. Dreyfus, R. D. Estes and R. E. Walkup, Appl. Phys. Lett., 1984, 44, 1155–1157.
83. G. Inoue and M. Suzuki, Chem. Phys. Lett., 1985, 122, 361–364.
84. D. Rayner, R. Steer, P. Hackett, C. Wilson and P. John, Chem. Phys. Lett., 1986, 123, 449–452.
85. J. Thoman and J. Steinfeld, Chem. Phys. Lett., 1986, 124, 35–38.
86. J. J. O'Brien and G. H. Atkinson, Chem. Phys. Lett., 1986, 130, 321–329.
87. J. W. Thoman, J. I. Steinfeld, R. I. McKay and A. E. W. Knight, J. Chem. Phys., 1987, 86, 5909–5917.
88. M. Fukushima, S. Mayama and K. Obi, J. Chem. Phys., 1992, 96, 44–52.
89. A. Kono, N. Koike, K. Okuda and T. Goto, Jpn. J. Appl. Phys., 1993, 32, L543–L546.
90. H. Nomura, K. Akimoto, A. Kono and T. Goto, J. Phys. D: Appl. Phys., 1995, 28, 1977–1982.
91. R. Escribano and A. Campargue, J. Chem. Phys., 1998, 108, 6249–6257.
92. S. N. Yurchenko, A. F. Al-Refaie and J. Tennyson, Astron. Astrophys., 2018, 614, A131.
93. S. N. Yurchenko, I. Szabo, E. Pyatenko and J. Tennyson, Mon. Not. R. Astron. Soc., 2018, 480, 3397–3411.
94. H. Y. Li, S. N. Yurchenko and J. Tennyson, Mon. Not. R. Astron. Soc., 2019, 486, 2351–2365.
95. A. Y. Adam, P. Jensen, A. Yachmenev and S. N. Yurchenko, J. Mol. Spectrosc., 2019, 362, 77–83.
96. A. N. Smirnov, V. G. Solomonik, S. N. Yurchenko and J. Tennyson, Phys. Chem. Chem. Phys., 2019, 21, 22794–22810.
97. K. A. Peterson, T. B. Adler and H.-J. Werner, J. Chem. Phys., 2008, 128, 084102.
98. S. Ten-No, Chem. Phys. Lett., 2004, 398, 56–61.
99. J. G. Hill, K. A. Peterson, G. Knizia and H.-J. Werner, J. Chem. Phys., 2009, 131, 194105.
100. K. E. Yousaf and K. A. Peterson, J. Chem. Phys., 2008, 129, 184108.
101. F. Weigend, Phys. Chem. Chem. Phys., 2002, 4, 4285–4291.
102. C. Hättig, Phys. Chem. Chem. Phys., 2005, 7, 59–66.
103. K. Tonokura, T. Murasaki and M. Koshi, J. Phys. Chem. B, 2002, 106, 555–563.
104. B. Beagley, A. R. Conrad, J. M. Freeman, J. J. Monaghan, B. G. Norton and G. C. Holywell, J. Mol. Struct., 1972, 11, 371–380.
105. J. Duncan, J. Harvie, D. McKean and S. Cradock, J. Mol. Struct., 1986, 145, 225–242.
106. M. D. Hanwell, D. E. Curtis, D. C. Lonie, T. Vandermeersch, E. Zurek and G. R. Hutchison, J. Cheminformatics, 2012, 4, 17.
107. D. L. Kokkin, T. Ma, T. Steimle and T. J. Sears, J. Chem. Phys., 2016, 144, 244304.
108. B. T. Sutcliffe and J. Tennyson, Int. J. Quantum Chem., 1991, 39, 183–196.
109. S. N. Yurchenko, T. M. Mellor, R. S. Freedman and J. Tennyson, Mon. Not. R. Astron. Soc., 2020, 496, 5282–5291.
110. P. Jensen, Mol. Phys., 2000, 98, 1253–1285.
111. A. Yachmenev and S. N. Yurchenko, J. Chem. Phys., 2015, 143, 014105.
112. S. N. Yurchenko, A. Yachmenev and R. I. Ovsyannikov, J. Chem. Theory Comput., 2017, 13, 4368–4381.
113. J. Tennyson and S. N. Yurchenko, Int. J. Quantum Chem, 2017, 117, 92–103.
114. P. R. Bunker and P. Jensen, Molecular Symmetry and Spectroscopy, NRC Research Press, Ottawa, 2nd edn, 1998.
115. B. V. Noumerov, Mon. Not. R. Astron. Soc., 1924, 84, 592–602.
116. J. W. Cooley, Math. Comp., 1961, 15, 363–374.
117. P. Jensen, M. Brumm, W. Kraemer and P. Bunker, J. Mol. Spectrosc., 1995, 171, 31–57.
118. S. N. Yurchenko, P. R. Bunker, W. P. Kraemer and P. Jensen, Can. J. Chem., 2004, 82, 694–708.
119. A. Campargue, D. Romanini and N. Sadeghi, J. Phys. D: Appl. Phys., 1998, 31, 1168–1175.

---

Footnote

† Electronic supplementary information (ESI) available. See DOI: 10.1039/d1cp00839k

---