A variationally computed room temperature line list for AsH$_{3}$

TL;DR

This paper presents a highly accurate variational calculation of the rotation-vibration energy levels and transition intensities of arsine (AsH₃), resulting in a comprehensive line list at room temperature that aligns well with experimental data.

Contribution

The study introduces a refined potential energy surface and a new dipole moment surface for AsH₃, enabling precise computation of energy levels and transition intensities up to 10,500 cm⁻¹ and J=30.

Findings

Energy levels match experimental data within 0.122 cm⁻¹ RMS error.

Calculated line intensities agree within 10% with experimental measurements.

Absorption cross-sections align well with laboratory spectra.

Abstract

Calculations are reported on the rotation-vibration energy levels of the arsine molecule with associated transition intensities. A potential energy surface (PES) obtained from \textit{ab initio} electronic structure calculations is refined to experimental data, and the resulting energy levels display sub-wavenumber accuracy for all reliably known $J=0$ term values under 6500 cm$^{-1}$. After a small empirical adjustment of the band centres, our calculated ($J=1-6$) rovibrational states reproduce 578 experimentally derived energies with a root-mean-square error of 0.122 cm$^{-1}$. Absolute line intensities are computed using the refined PES and a new dipole moment surface (DMS) for transitions between states with energies up to 10~500 cm$^{-1}$ and rotational quantum number $J=30$. The computed DMS reproduces experimental line intensities to within 10\% uncertainty for the $\nu_1$ and $\nu_3$ bands. Furthermore, our calculated absorption cross-sections display good agreement with the main absorption features recorded in Pacific Northwest National Laboratory (PNNL) for the complete range of $600-6500$ cm$^{-1}$.

Tables

Table 1. Table 1 : Equilibrium energies calculated at the CCSD(T) level of theory using different basis sets and Hamiltonians.

Basis	Energy/$E_h$
AVQZ	-2236.17795527
AVQZ-DK	-2261.02359376
AVQZ-PP	-333.14700414
AV5Z	-2236.18098311
AV5Z-DK	-2261.02786843
AV5Z-PP	-333.14983002

Table 2. Table 2 : Differences between experimentally derived band centres and our calculated values computed using all-electron DKH and pseudopotential-F12 based PESs. All numerical values are energies given in units of cm -1 .

Band	Symmetry	Obs.	VQZ-PP-F12	AVQZ-DK	Obs$-$CalcPP-F12	Obs$-$CalcDK
${\nu }_2$	A1	906.752	904.812	905.058	1.940	1.694
${\nu }_4$	E	999.225	994.460	994.132	4.765	5.093
$2{\nu }_2$	A1	1806.149	1802.443	1802.451	3.706	3.698
${\nu }_2+{\nu }_4$	E	1904.115	1897.551	1897.465	6.564	6.650
$2{\nu }_4^{l=0}$	A1	1990.998	1982.116	1981.574	8.882	9.424
$2{\nu }_4^{l=2}$	E	2003.483	1988.246	1987.651	15.237	15.832
${\nu }_1$	A1	2115.164	2108.659	2105.000	6.505	10.164
${\nu }_3$	E	2126.432	2116.469	2112.542	9.963	13.890
${\nu }_1+{\nu }_2$	A1	3013a	3006.718	3002.875	6.3	10.1
${\nu }_1+{\nu }_4$	E	3102a	3089.255	3084.866	12.7	17.1
$2{\nu }_1$	A1	4166.772	4151.833	4143.187	14.939	23.585
${\nu }_1+{\nu }_3$	E	4167.935	4152.229	4143.527	15.706	24.408
$2{\nu }_3^{l=0}$	A1	4237.700	4222.006	4214.312	15.694	23.388
$2{\nu }_3^{l=2}$	E	4247.720	4229.805	4221.816	17.915	25.904
$2{\nu }_1+{\nu }_2$	A1	5057a	5041.541	5030.916	15.5	26.1
${\nu }_1+{\nu }_2+{\nu }_3$	E	5057a	5041.191	5030.920	15.8	26.1
$2{\nu }_1+{\nu }_4$	E	5128a	5111.286	5100.477	16.7	27.5
$2{\nu }_3^0+{\nu }_2$	A1	5128a	5113.615	5104.249	14.4	23.8
${\nu }_1+{\nu }_3+{\nu }_4$	E	5158a	5137.282	5127.176	20.7	30.8
${\nu }_1+{\nu }_3+{\nu }_4$	A1	5158a	5137.555	5127.471	20.4	30.6
$3{\nu }_1$	A1	6136.340	6116.822	6101.231	19.518	35.109
$2{\nu }_1+{\nu }_3$	E	6136.330	6116.793	6101.192	19.537	35.138
${\nu }_1+2{\nu }_3^{l=0}$	A1	6275.830	6257.116	6243.540	18.714	32.290
${\nu }_1+2{\nu }_3^{l=2}$	E	6282.350	6261.282	6247.600	21.068	34.750
$3{\nu }_3^{l=1}$	E	6294.710	6270.037	6256.059	24.673	38.651
$3{\nu }_3^{l=3}$	A1	6365.950	6340.980	6327.902	24.970	38.048

Table 3. Table 3 : Experimental and predicted structural constants of 75 AsH 3 .

	AsH3-CYT18	VQZ-PP-F12	AVQZ-DK	AV5Z-DK	Exp 65
$r_{\mathrm{eq}}$ /Å	1.511394	1.520269	1.521481	1.520432	1.511060
${\alpha }_{\mathrm{eq}}$ /∘	92.04025	92.21595	92.17049	92.18705	92.0690
$r_{\mathrm{SP}}$ /Å	1.4688	1.4663	1.4670
$\Delta E$(barrier) / cm-1	14495.	14187.	14171.

Table 4. Table 4 : Differences between calculated rotational energies, in cm -1 , and the hyperfine resolved values of 62 which we averaged using the spin statistical weights.

J	K	Sym	Obs	Obs$-$Calcref
1	0	A2	7.503018	-0.000368
1	1	E	7.249824	-0.000381
2	1	E	22.253842	-0.001287
2	2	E	21.494930	-0.000807
3	0	A2	45.005427	-0.003161
3	1	E	44.753718	-0.002935
3	2	E	43.997254	-0.002248
3	3	A1	42.732027	-0.001078
3	3	A2	42.732025	-0.001077
4	1	E	74.742676	-0.005085
4	2	E	73.989211	-0.004375
4	3	A1	72.729005	-0.003149
4	3	A2	72.728988	-0.003184
4	4	E	70.955354	-0.001422
5	0	A2	112.460919	-0.007935
5	1	E	112.211432	-0.007705
5	2	E	111.461647	-0.007010
5	3	A1	110.207627	-0.005823
5	3	A2	110.207557	-0.005824
5	4	E	108.442621	-0.004109
5	5	E	106.157387	-0.001816
6	1	E	157.148421	-0.010773
6	2	E	156.403009	-0.010103
6	3	A1	155.156374	-0.008751
6	3	A2	155.156164	-0.009170
6	4	E	153.401598	-0.007307
6	5	E	151.129706	-0.005092

Table 5. Table 5 : Agreement between our calculated energy levels and those derived from experiment. All calculations used our refined PES, AsH 3 -CYT18. J = 0 𝐽 0 J=0 comparisons are before employing the EBSC, and J = 1 − 6 𝐽 1 6 J=1-6 comparisons are afterwards. Energy units are cm -1 .

Band	Symmetry	$J=0$				$J=1-6$
			Obs.	Calc.	Obs.$-$Calc.		${\sigma }_{\mathrm{rms}}^{\mathrm{ebsc}}$
${\nu }_2$	A1		906.752	906.109	0.643		$-$
${\nu }_4$	E		999.225	998.833	0.393		$-$
$2{\nu }_2$	A1		1806.149	1806.161	-0.012		0.048
${\nu }_2+{\nu }_4$	E		1904.115	1904.046	0.069		0.131
$2{\nu }_4^{l=0}$	A1		1990.998	1990.293	0.705		0.262
$2{\nu }_4^{l=2}$	E		2003.483	1997.315	6.168		0.207
${\nu }_1$	A1		2115.164	2114.938	0.227		0.027
${\nu }_3$	E		2126.432	2126.102	0.330		0.068
${\nu }_1+{\nu }_2$	A1		3013a	3016.531	-3.5		$-$
${\nu }_1+{\nu }_4$	E		3102a	3100.438	2.4		$-$
$2{\nu }_1$	A1		4166.772	4166.694	0.078		0.067
${\nu }_1+{\nu }_3$	E		4167.935	4167.877	0.058		0.059
$2{\nu }_3^{l=0}$	A1		4237.700	4237.407	0.293		0.046
$2{\nu }_3^{l=2}$	E		4247.720	4247.842	-0.122		0.241
$2{\nu }_1+{\nu }_2$	A1		5057a	5040.690	16.3		$-$
${\nu }_1+{\nu }_2+{\nu }_3$	E		5057a	5040.799	16.2		$-$
$2{\nu }_1+{\nu }_4$	E		5128a	5129.956	-2.0		$-$
$2{\nu }_3^0+{\nu }_2$	A1		5128a	5131.122	-3.0		$-$
${\nu }_1+{\nu }_3+{\nu }_4$	E		5158a	5155.741	2.3		$-$
${\nu }_1+{\nu }_3+{\nu }_4$	A1		5158a	5156.434	1.6		$-$
$3{\nu }_1$	A1		6136.340	6136.846	-0.506		$-$
$2{\nu }_1+{\nu }_3$	E		6136.330	6136.859	-0.529		$-$
${\nu }_1+2{\nu }_3^{l=0}$	A1		6275.830	6275.814	0.017		0.051
${\nu }_1+2{\nu }_3^{l=2}$	E		6282.350	6282.414	-0.064		0.049
$3{\nu }_3^{l=1}$	E		6294.710	6294.695	0.015		0.050
$3{\nu }_3^{l=3}$	A1		6365.950	6365.759	0.191		0.030

Table 6. Table 6 : Comparison of observed and calculated band intensities. Column 1 refers to the local mode quantum numbers assigned by TROVE, where sym is the total symmetry. The units of intensity are 10 -18 cm / − 1 ( molec {}^{-1}/(\mathrm{molec} cm -2 ). The value under / is the total intensity of the bands with the same quantum numbers n 1 ​ n 2 ​ n 3 subscript 𝑛 1 subscript 𝑛 2 subscript 𝑛 3 n_{1}n_{2}n_{3}

($n_1,n_2,n_3$;sym)	band	bend centre	$I_{\mathrm{obs}}$71	$I_{\mathrm{calc}}$71	$I_{\mathrm{calc}}$ (this work)
(100;$A_1$)	${\nu }_1$	2115.164	11.4/44.1	10.7/40.4	11.2/44.9
(100;$E$)	${\nu }_3$	2126.432	29.7	32.7	33.7
(200;$A_1$)	$2{\nu }_1$	4166.772	/0.618	0.157/0.427	0.225/0.722
(200;$E$)	${\nu }_1+{\nu }_3$	4167.935	$-$	0.270	0.497
(110;$A_1$)	$2{\nu }_3^{l=0}$	4237.700	$-$	0.0117/0.0123	0.0143/0.0163
(110;$E$)	$2{\nu }_3^{l=2}$	4247.720	$-$	0.000671	0.00201
(300;$A_1$)	$3{\nu }_1$	6136.340	/0.00989	0.00456/0.00656	0.00337/0.00548
(300;$E$)	$2{\nu }_1+{\nu }_3$	6136.330	$-$	0.00200	0.00211
(210;$A_1$)	${\nu }_1+2{\nu }_3^{l=0}$	6275.830	/0.00275	0.00112/0.00182	0.000734/0.00116
(210;$E$)	${\nu }_1+2{\nu }_3^{l=2}$	6282.350	$-$	0.000104	0.0000741
(210;$E$)	$3{\nu }_3^{l=1}$	6294.710	$-$	0.000596	0.000356
(111;$A_1$)	$3{\nu }_3^{l=3}$	6365.950	$-$	0.0000650	0.0000934

Table 7. Table 7 : Comparison of calculated and observed 15 line positions and intensities belonging to the ν 1 subscript 𝜈 1 \nu_{1} and ν 3 subscript 𝜈 3 \nu_{3} bands.

$J$’	$K^{\prime }$	Sym′	$J^{\prime \prime }$	$K^{\prime \prime }$	Sym′′	band	${\nu }_{\mathrm{obs}}$15	$I_{\mathrm{obs}}$15	${\nu }_{\mathrm{calc}}$	$I_{\mathrm{calc}}$	%$\left|\frac{I_{\mathrm{obs}}-I_{\mathrm{calc}}}{I_{\mathrm{obs}}}\right|$
9	6	E	10	7	E	${\nu }_3$	2051.894	4.799$\times$10-20	2052.082	4.992$\times$10-20	4.03
9	7	E	10	8	E	${\nu }_3$	2052.548	5.755$\times$10-20	2052.767	6.079$\times$10-20	5.63
7	1	A2	8	0	A1	${\nu }_3$	2064.460	4.396$\times$10-20	2064.468	4.773$\times$10-20	8.57
7	6	E	8	7	E	${\nu }_3$	2067.961	1.038$\times$10-19	2068.139	1.110$\times$10-19	6.95
4	4	E	5	5	E	${\nu }_3$	2090.433	1.572$\times$10-19	2090.542	1.699$\times$10-19	8.04
4	1	A1	5	0	A2	${\nu }_3$	2088.098	5.876$\times$10-20	2088.087	6.423$\times$10-20	9.31
5	4	E	6	5	E	${\nu }_3$	2082.601	1.250$\times$10-19	2082.714	1.320$\times$10-19	5.64
3	3	E	4	4	E	${\nu }_3$	2097.659	1.462$\times$10-19	2097.738	1.571$\times$10-19	7.49
2	2	E	1	1	E	${\nu }_3$	2140.716	8.659$\times$10-20	2140.678	9.453$\times$10-20	9.17
2	1	A1	1	0	A2	${\nu }_3$	2141.069	8.949$\times$10-20	2141.053	9.612$\times$10-20	7.41
6	6	E	5	5	E	${\nu }_3$	2168.331	1.965$\times$10-19	2168.331	2.093$\times$10-19	6.48
8	5	E	7	5	E	${\nu }_1$	2172.196	3.114$\times$10-20	2172.262	3.301$\times$10-20	5.99
8	7	E	9	7	E	${\nu }_1$	2045.190	2.459$\times$10-20	2045.261	2.326$\times$10-20	5.40
8	8	E	9	8	E	${\nu }_1$	2045.319	1.310$\times$10-20	2045.397	1.328$\times$10-20	1.36
10	7	E	9	7	E	${\nu }_1$	2185.605	2.026$\times$10-20	2185.716	2.114$\times$10-20	4.35

Equations

$$V({\xi }_1,{\xi }_2,{\xi }_3,$$

$$+\sum _{i\le j}{F}_{ij}\sin (\bar{\rho }){\xi }_i{\xi }_j+...$$

$$+\sum _{i\le j\le k\le l\le m\le n}{F}_{ijklmn}\sin (\bar{\rho }){\xi }_i{\xi }_j{\xi }_k{\xi }_l{\xi }_m{\xi }_n$$

$${\xi }_k$$

$${\xi }_4$$

$${\xi }_5$$

$$\sin \bar{\rho }$$

$$F_{ij...}\sin (\bar{\rho })=\sum _{s=0}^N{f}_{ij...}^{(s)}(\sin ({\rho }_{\mathrm{eq}})-\sin (\bar{\rho }){)}^s$$

$$w_i=\frac{2}{1+e^{2\times 10^{-4}\times E_i}}$$

$$|\nu ,J,K,m,{\tau }_{\mathrm{rot}}⟩=[|J,K,m,{\tau }_{\mathrm{rot}}⟩|n_1⟩|n_2⟩|n_3⟩|n_4⟩|n_5⟩|n_6⟩{]}^{{\Gamma }_{\mathrm{ir}}}$$

$$\mathrm{P}=2(n_1+n_2+n_3)+n_4+n_5+n_6$$

$$S=$$

$$+k\sum _m{w}_m\{E_m^{\text{ab}}-E_m^{ref}(f_{jk...}+\Delta f_{jk...}){\}}^2$$

$$H=T+V+\Delta V$$

$${\bar{\mu }}_X=\frac{\partial E}{\partial F_X}=\frac{E(+\Delta F_X)-E(-\Delta F_X)}{2\Delta F_X},$$

$${\bar{\mu }}_Y=\frac{\partial E}{\partial F_Y},$$

$${\bar{\mu }}_Z=\frac{\partial E}{\partial F_Z}.$$

$${\overline{\mu }}_{\Gamma }({\chi }_1,{\chi }_2,{\chi }_3,{\chi }_{4a},{\chi }_{4b};\rho )=$$

$$+\sum _{i\le j\le k}{\mu }_{ijk}^{\Gamma }(\rho ){\chi }_i{\chi }_j{\chi }_k$$

$${\chi }_k$$

$${\chi }_4$$

$${\chi }_5$$

$${\mu }_{ij...}^{\Gamma }(\rho )=\sum _{s=0}^N{\mu }_{ij...}^{\Gamma (s)}(\sin ({\rho }_e)-\sin (\rho ){)}^s,$$

$$I(f\leftarrow i)=$$

$$\times \sum _{{\Phi }_{int}^{\prime }{\Phi }_{int}^{\prime \prime }}\sum _{A=X,Y,Z}|⟨{\Phi }_{int}^{\prime }|{\mu }_A|{\Phi }_{int}^{\prime \prime }⟩{|}^2$$

$$Q=\sum _w{g}_w\exp (-E_w/kT)$$

Full text

\makeFNbottom

[TABLE]

††footnotetext: a*  Department of Physics & Astronomy, University College London, London WC1E 6BT, UK ††footnotetext: b  Servomex Ltd., Millbrook Industrial Estate, Crowborough TN6 3FB, UK *††footnotetext: † Electronic Supplementary Information (ESI) available: including the final line list, partition function and other spectroscopic data is available from the ExoMol website, www.exomol.com. See DOI: 10.1039/b000000x/

1 Introduction

Arsine (AsH3) is a highly poisonous gas1 which is the direct analogue molecular structure of ammmonia (NH3) and phospine (PH3). Like these two gases it has been detected the atmospheres of the gas giant planets Jupiter 2, 3 and Saturn 4. It may therefore be expected to be also present in the atmospheres of gas giant exoplanets.

Arsene is also important for industrial applications as high purity arsine is widely used in the semiconductor manufacturing industry, for example, in processing GaAs surfaces 5, 6, 7. Given its highly poisonous nature, with an exposure limit value of 50 ppb mole-concentration8, the detection of AsH3 escape at such levels is an important safety requirement in this industry9. It is also monitored in the polymer industry as trace level arsine impurity in ethylene and propylene monomer feedstock gases may contaminate the catalysts, resulting in reduced quality and yield of the polymer products 10.

Arsine is also a trace atmospheric pollutant due to emissions from various industrial processes, such as power generation and smelting 11, 12. Routine methods for arsine measurement in industry include gas chromatography, electrochemical sensors, colorimetric sensors, and Fourier-transform infrared spectroscopy. However, the development of high resolution laser spectroscopy based measurements is a growing area13, 14 for which detailed line lists are required to model the high resolution absorption spectra.

While there have been a number of studies of the infrared and microwave spectrum of arsine, there is no comprehensive line list for the system and there is a lack of information on the intensity of many bands. The situation for absolute line intensities is particularly dire, with existing data confined solely to the measurements reported by Dana et al. 15. Previous attempts to model the global vibrational structure 16, 17, and rovibrational sub-structures18, 19, 20, 21, 22, 23, 24, have focussed predominantly on Effective-Hamiltonians, which have limited predictive capability outside the fitted data. In addition, an ab initio potential energy surface (PES) for AsH3 was reported in 1995 by Breiding and Thiel in the form of the cubic anharmonic force field 25 using relativistic effective core potentials (ECPs).

Considering the unsuitability of the current state of AsH3 data for either exoplanet modelling, which necessitates completeness, or industrial monitoring, which necessitates accuracy, we decided to construct a comprehensive line list for arsine which could be used for the applications mentioned above. Our approach to constructing linelists, as exploited in the ExoMol project 26, 27, uses potential energy surfaces which have been refined using spectroscopic data but ab initio dipole moment surfaces, which have been shown to give highly accurate predicted transition intensives 28, 29, 30. The approach has already been used to compute line lists for the systems NH3 31, 32, 33, PH3 34, 35, SbH3 36, SO3 37 and and a linelist for PF3 is on the way 38. These studies were based on either a D3h symmetry 31, 32, 33, 39, in which the tunneling of the so-called umbrella motion was explicitly considered, or C3v symmetry 34, 35 in which the tunneling is neglected. Given the high barrier expected for AsH3, the reduced symmetry C3v approach is employed here.

2 Potential energy surface

2.1 Electronic structure calculations

Accurate modelling of heavy elements in quantum chemistry is made particularly challenging by increased relativistic effects, core-core electron correlation and core-valence electron correlation compared to lighter elements. Recovery of the correlation energy has benefited from the recent development of explicitly correlated methods (F12/R12) 40, 41, which rapidly converge electronic energies towards the complete basis set (CBS) limit with increasing cardinal number $n$, but must be used in conjunction with suitably optimised basis sets for full effectiveness. To account for scalar relativistic effects, effective core potentials (ECPs) are a computationally inexpensive option. More rigorous treatment is possible with the Douglas-Kroll-Hess (DKH) Hamiltonian, however, so far no satisfactory way has been found of incorporating the DKH Hamiltonian into the F12/R12 framework 42. Although ECPs face the same theoretical hurdles, namely non-commutability with the F12 correlation function, alternative treatments have been found to work well 43, 44. For heavy elements where both F12-pseudopotential and standard all-electron DKH based approaches are possible, such as arsenic, the benefits of F12 must therefore be weighed against the penalty of introducing an additional scalar-relativistic approximation. Peterson 45 showed that complete basis set (CBS) extrapolated CCSD(T)/aug-cc-pwCVnZ-PP pseudopotential calculations performed almost identically to their DKH all-electron counterparts in a series of molecular benchmark calculations for post-3d main group elements, including the As2, AsF, AsCl and AsN molecules. They go on to develop a new family of F12-specific cc-pVnZ-PP-F12 basis sets to be used at the CCSD(T)-F12 level, which yield accuracy comparable to the 2-3 times larger aug-cc-pwCV(n+2)Z-PP basis sets used at standard CCSD(T) level 45, 46. Their pseudopotential-F12 optimised approach is the one followed in this work.

All electronic structure calculations were performed using MOLPRO 47 and employed the explicitly correlated coupled cluster method CCSD(T)-F12b 40, 41 with implicit treatment of scalar-relativistic effects via replacement of 10 core electrons with a pseudopotential (PP). Calculations were carried out in the frozen core approximation and utilized the correlation consistent quadruple-zeta, PP-F12 optimised basis set of Hill et al. 46 (cc-pVQZ-PP-F12) to represent the As electronic wavefunction, and cc-pVQZ-F12 basis sets for the H atoms. Density fitting (DF) for the 2-electron (MP2FIT) and exchange term (JKFIT) integrals employed the cc-pVTZ-PP-F12/MP2Fit and def2-QZVPP/JKFIT basis sets respectively, and for the resolution of the identity of the many-electron F12 integrals (OPTRI) we used the VTZ-PP-F12/OPTRI basis set. For the geminal exponent $\gamma$, a value of 1.4 a${}_{0}^{-1}$ was used as recommended by Hill et al.. All calculations were performed on the ground electronic state, which is sufficiently uncoupled from higher electronic excitations that both adiabatic and non-adiabatic effects are expected to be very small 48.

To initially qualify the importance of including scalar relativistic effects in our calculations, relativistic corrections $\Delta E_{X}$ along 1-dimensional cuts of the potential energy surface PES were calculated (shown in Fig. 1, along with cuts through the AVQZ surface for reference). This was done by by first shifting the potential energy curves by their respective energies at equilibrium, listed in table 1, where we note the relatively small absolute energies of the ECP based calculations owing to their implicit treatment of 10 core electrons. We then have $\Delta E_{X}=E_{X_{1}}-E_{X_{2}}$ where $E_{X_{1}}=$AV$n$Z-PP/AV$n$Z-DK and $E_{X_{2}}=$AV$n$Z for the pseudopotential/all-electron calculations. Here, and in all subsequent DKH calculations, the DKH Hamiltonian has been expanded to 8th–order (DKH8) using optimal unitary parametrisation. In Fig. 1 only quadruple-zeta ($n=4$) results are presented as they were seen to differ by no more than 3 cm*-1* from the respective 5-zeta ($n=5$) curves. Clearly the inclusion of scalar relativistic effects are important, and both approximations have a similar effect on the total energy. However, the pseudopotential approximation tends to raise the energy at stretched geometries and lower the energy at contracted geometries, relative to the all-electron calculations. It is difficult to asses the effect of this difference within the Hill et al. 46 regime for a full-dimensional PES. We therefore opted to generate a second 6D surface at the AVQZ-DKH8 level of theory (henceforth denoted AVQZ-DK), to provide a benchmark for our VQZ-PP-F12 based ab initio nuclear motion calculations. The results of these are presented in section 3.

It is known that beyond 4th–order expansion the DKH Hamiltonian depends slightly on the chosen paramterization of the unitary transformations applied49. Thus, additional DKH4 calculations were performed on a subset of the DKH8 PES geometries. Between 4th and 8th–order the resulting electronic energies were seen to differ by less than 0.1 cm*-1* above their equilibrium values. For the purpose of benchmarking the ECP based nuclear motion calculations, therefore, this dependency is not expected to be significant.

2.2 Nuclear geometry grid

Our grid of nuclear geometries was built by combining 1D$-$6D sub-grids. Our 1D grid consisted of a cut along the $r_{1}=r_{2}=r_{3}$ stretch with $\alpha_{1}=\alpha_{2}=\alpha_{3}=\alpha_{eq}$, and a cut along the $\alpha_{1}=\alpha_{2}=\alpha_{3}$ bend with $r_{1}=r_{2}=r_{3}=r_{\rm eq}$. Each additional degree of freedom was then added by allowing either an additional As-H bond length or H-As-H bond angle to vary. Because this method causes the number of points to grow so rapidly, the range and intervals of $r$ and $\alpha$ were reduced with each increasing degree of freedom using the 1D cuts as a guide. This also helped to limit the range of electronic energies generated, as large distortions in geometry can lead to unnecessarily high values of energy that are not needed in the fit.

In order to ensure each grid point was fully unique we applied the C3v molecular symmetry group transformations prior to computing the electronic energy. If two grid points were transformed into one another, then one was discarded. Finally, any energetically sparse regions were filled by generating additional geometries that were estimated to fall within our desired range. The 1D cuts provided an initial guide to the electronic energy, then intermediate versions of our PES were used to more accurately choose geometries. Our final grid consisted of 39 873 nuclear geometries within the range $1.10\leq r_{i}\leq 3.74$ Å and $37^{\circ}\leq\alpha_{i}\leq 130^{\circ}$, with electronic energies extending to 27 000 cm*-1*, although $\sim$ 38 000 of these were below 10,000 cm*-1* ($1.25\leq r_{i}\leq 1.9$ Å and $60^{\circ}\leq\alpha_{i}\leq 126^{\circ}$). The additional points in the 10 000$-$27 000 cm*-1* range predominantly belonged to the As-H dissociative stretch, which is where holes commonly appear if the function is not suitably constrained at high energy. Grid points for our AVQZ-DK reference PES were chosen by randomly sampling 16 396 equally energetically distributed points from our VQZ-PP-F12 grid, which spanned the bond lengths $1.3\leq r_{i}\leq 1.8$ Å, bond angles $65^{\circ}\leq\alpha_{i}\leq 130^{\circ}$ and energies below 15 000 cm*-1*.

Each grid point computed at the CCSD(T)-F12b/cc-pVQZ-PP-F12 level took approximately 10-15 minutes to compute on UCL’s Legion computer cluster. This was increased to 20-30 mins for the DKH Hamiltonian-based calculations, owing to the increased computational demand of explicitly treating the 10 core electrons.

2.3 Analytic representation

To represent the PES analytically we used the same functional form as for NH3 and PH3 50. The potential is represented as a polynomial expansion

[TABLE]

in terms of the internal coordinates

[TABLE]

In Eq. (1)

[TABLE]

and $r_{k}$ is the As-Hk bond length, $\alpha_{j}$ is the $j^{th}$ H-As-H bond angle (opposite to the $j$th bond), $r_{\rm eq}$ is the equilibrium value of $r_{k}$, $a$ is a molecular parameter, and $\rho_{\rm eq}$ is the equilibrium value of the umbrella mode $\bar{\rho}$. $V_{0}$ represents the pure inversion potential and $f_{ij...}^{(s)}$ are the fitting parameters. Points were given energy ($E_{i}$) dependant weights ($w_{i}$)

[TABLE]

as used by Polyansky et al. 51. We could usefully fit terms in the potential up to 5th order resulting in a root-mean-square (RMS) deviation of 0.7 cm*-1* for the 39 678 nuclear geometries. For our all-electron reference PES, the weighted RMS error increased to 1.2 cm*-1*, most likely due to the proportionally fewer points close to equilibrium. However it should be noted that adding more points to the fit had little effect on the computed vibrational term values reported in section 3.

The final set of expansion parameters for our ab initio PES is included in the supplementary material, along with our grid of nuclear geometries and an Fortran 90 routine to evaluate the analytic expression.

3 Nuclear motion calculations

To calculate rovibrational energy levels we used the variational nuclear motion program TROVE. The general methodology of TROVE is well documented 52, 53, with its specific application to XY3-type molecules detailed in Ref.  31, and so only a brief description is provided here.

Rovibrational basis functions are constructed as symmetrised linear combinations of 1D primitive-basis-function products

[TABLE]

where the 1D stretching functions ($|n_{1}\rangle$,$|n_{2}\rangle$,$|n_{3}\rangle$) and bending functions ($|n_{4}\rangle$,$|n_{5}\rangle$,$|n_{6}\rangle$) are obtained by solving the corresponding one-dimensional Schrödinger equations using the Numerov-Cooley approach54, 55 for the stretches, and 1D harmonic oscillator eigenfunctions for the bends. The rotational basis set is built from rigid-rotor functions. In the above equation, $\Gamma_{\rm ir}$ represents one of the irreducible representations of C3v spanned by $|\nu,J,K,m,\tau_{\rm rot}\rangle$. For details of the symmetrisation procedure in TROVE the reader is directed to Ref. 53. A multi-step contraction scheme was employed to limit the vibrational, then rovibrational basis set size. This is outlined in the following paragraphs and in section 6.2.

Owing to the structural similarities between AsH3 and other XY3-type molecules which have been investigated in the past, variational calculations could be performed with relative ease once a PES and DMS had been constructed, and required only a few molecule specific parameters to be defined: atomic mass, molecular symmetry group, Z-matrix, equilibrium parameters, and the definition of our 1-D grids upon which the wavefunctions are evaluated. In fact, this is true for any 3-5 atom molecule, providing the molecular symmetry (MS) group has been programmed into TROVE, this includes a recent extension to treat linear molecules 56. Currently TROVE allows for the following MS groups: C2v, C3v, D3h, D2d, C2h, Td, D2h, as well as Dnh with $n>1$ 57 (see also Ref.  53, where the TROVE treatment of the symmetry groups is described). Here we use C3v symmetry which means that the inversion mode is not fully represented and any tunneling splitting is neglected. We note that study of deuterated arsine, AsH2D could also be conducted using the same Born-Oppenheimer PES but would require a different symmetrisation procedure in TROVE. Such a study is possible but the lower symmetry and denser vibrational spectrum associated with the D atom would make such a calculation computationally more expensive.

As arsenic has only one stable isotope, 75As, all reported nuclear motion calculations were performed for 75AsH3.

Within the limitations of our PES, the accuracy of our variational calculation is determined predominantly by i) the size of our nuclear-motion basis set; and ii) our Taylor-type expansion of the kinetic energy operator $\hat{T}$ and our re-expansion of the potential function $V$ in terms of linearised coordinates. The former, we choose to restrict via a polyad number P. The polyad number is an integer that represents the total quanta of vibrational excitations in terms of the lowest energy fundamental. For AsH3

[TABLE]

where $n_{1}+n_{2}+n_{3}$ is the total number of stretching quanta and $n_{4}+n_{5}+n_{6}$ is the total number of bending quanta, corresponding to the primitive functions $|n_{i}\rangle$ ($i=1,...,6$) in Eq.(8). For our comparison of the VQZ-PP-F12 and AVQZ-DK based ab initio PESs we chose to include in our variational calculations all vibrational states with P $\leq$ P${}_{\mathrm{max}}=$14 as used previously for NH3 and PH3 31, 34. This resulted in our vibrational eigenvalues converged to within 0.1 cm*-1* below 6000 cm*-1* for the stretches, and as much as 3 cm*-1* for the bending overtones. Our $\hat{T}$ and $V$ expansions we take to 6th and 8th order respectively. Increasing these to 8th and 10th order changed the vibrational term values reported throughout this work by $<0.1$ cm*-1* for the stretches, and $<0.7$ cm*-1* for the bends. For highly excited bending overtones, such as the 5 and 6–quanta bends, the convergence error due to our $\hat{T}$ and $V$ expansions may be several wavenumbers.

Table 2 shows the 26 lowest-lying experimentally derived vibrational states compared to our calculations. Term values known to sub-wavenumber accuracy are taken from Sanzharov et al 16; the remaining eight bands are from the work by Halonen et al. 58 and have an estimated 2 cm*-1* uncertainty, although this may be larger for the 5050$-$5200 cm*-1* bands 59. Using the VQZ-PP-F12 and AVQZ-DK PESs the experimentally derived values of the four fundamentals are reproduced to within 10 cm*-1* and 14 cm*-1* respectively. Whilst far from the accuracy achieved in previous studies of NH3 and PH3, our results are comparable to the achievements of Nikitin et al. in their recent ab initio study of GeH4 60, and we deem it reasonable considering the greater contribution of relativistic effects, core-core electron correlation and core-valence electron correlation associated with heavier atoms. For the overtones and combination bands the quality of our ab initio predictions steadily decreases in proportion to the error on the fundamentals, except for the $2\nu_{4}^{l=0}$ band which is independently examined in section 6.2. Most importantly, the VQZ-PP-F12 surface consistently and significantly outperforms the VQZ-DK surface. Given the factor of 2 reduction in computational time, this highlights the value of the work by Hill, Peterson and co-authors 46, 45.

4 Refinement

In order to achieve so-called ‘spectroscopic’ accuracy in our variational nuclear motion calculations it is common practice to empirically refine the chosen ab initio PES to experimental data. In this case our chosen starting point PES is the VQZ-PP-F12 surface. A theoretical description of the refinement procedure has been previously reported for the case of NH3 50, and our method for AsH3 is similar. Namely, small corrections to the parameters $f_{jk...}$ in equation (1) are made, so as to minimise the sum of squared residuals 61

[TABLE]

The energies $E_{n}^{calc}$ are found by diagonalising the matrix representation of the Hamiltonian

[TABLE]

where $\Delta V$ has the same form as Eq. (1), except $f_{jk...}$ are replaced by the adjustable parameters $\Delta f_{jk...}$. In the above equation $w_{i}$ are the weights applied, $E_{n}^{obs}$ are the experimentally derived energies, and $E_{m}^{\textit{ab}}$ and $E_{m}^{\textit{ref}}$ are the energies of the ab initio and refined PESs when evaluated on our grid of nuclear geometries. The second term in Eq. (10) ensures our refined potential retains the general shape of the ab initio surface, how strongly we force it to do so is controlled by the constant $k$. In order to find the set of parameters $\Delta f_{jk...}$ for which the above function is minimised, we employ an iterative least-squares fitting algorithm, which requires only the energies and their derivatives with respect to the adjustable parameters.

Because As is heavier than N or P, the rotational energies of AsH3 are more closely spaced than those of NH3 and PH3, and so more highly populated at room temperature. Particular attention was therefore paid to optimising the equilibrium bond lengths and bond angles. This optimisation was performed prior to the refinement by using the hyperfine resolved rotational energies of Tarrago et al. 62, which we averaged using the spin-statistical weights ({$A_{1}$,$A_{2}$,$E$}$=${16,16,16}), and a Newton-Gauss style procedure with a step size of $\pm 0.002$ Å and $\pm 0.002$ rad. Although TROVE is capable of computing quadrupole-hyperfine effects 63, requiring only a quadrupole moment surface and electric field gradient tensor in addition to the PES and DMS, the resulting splittings are small (roughly a few MHz) and so not considered here.

For the full nonlinear least squares refinement we allowed for corrections to harmonic and certain cubic terms in our PES, and used 322 experimentally derived energies with $J\leq 6$ compiled from Refs. 62, 18, 19, 20, 21, 22. These sampled the following vibrational bands: the fundamentals $\nu_{1}$, $\nu_{2}$, $\nu_{3}$, $\nu_{4}$; overtones $2\nu_{1}$, $2\nu_{2}$, $2\nu_{3}$, $2\nu_{4}$, $3\nu_{1}$, $3\nu_{3}$; and combination bands $\nu_{1}+\nu_{3}$, $\nu_{2}+\nu_{4}$, $2\nu_{1}+\nu_{3}$, $\nu_{1}+2\nu_{3}$. Because we could find no rotationally excited states belonging to the $\nu_{2}$ and $\nu_{4}$ bands in the literature, only their band centres were included in the refinement.The vibrational band centres measured by Halonen et al. 58 were not included due to the large estimated uncertainty. The complete list of experimental energies included in the refinement, along with their assignments, are included in the supplementary material.

Weights of $w_{n}=0.1$ were distributed to all experimentally derived rovibrational states except for purely rotational states, which were given weights of 1000.0, and the $2\nu_{4}^{l=0}$ band of Yang et al.20, for which we struggled to match experimental energies to our calculated energies owing to conflicting quantum labels, and so gave a weight of 0.0. These were adjusted on-the-fly using Watson’s robust fitting scheme 64. A scaling factor of $k=1\times 10^{-4}$ was initially applied to the 39 678 ab initio points $E_{m}^{\textit{ab}}$ included in the fit. As the refinement progressed this was incrementally decreased to $1\times 10^{-6}$ so as to reduce the relative contribution of the ab initio data. Care was taken throughout to ensure the refined PES did not deviate substantially from the ab initio surface, and we note that for all ab initio grid points, the energy difference between the refined and geometry optimised ab initio PES’s is less than 10$\%$ that of the ab initio PES above its zero-point energy (ZPE).

Our final fitted PES is called AsH3-CYT18 below; it is presented as subroutine in the supplementary material.

5 Dipole moment surface

5.1 Electronic structure calculations

The electric dipole moment is equal to the first derivative of the electronic energy with respect to an external electric field. This can be approximated using a numerical finite-difference procedure whereby the dipole moment $\overline{\mu}$ is related to the electronic energy in the presence of a weak uniform electric field $\Delta F$ acting along the space fixed $XYZ$ axis

[TABLE]

[TABLE]

[TABLE]

Thus to evaluate the dipole moment at a given nuclear geometry requires seven electronic structure calculations, the seventh being an initial zero-field calculation. A field strength of 0.002 a.u. was deemed sufficiently small to accurately approximate the first derivative without approaching numerical noise 30. As with our PES, electronic structure calculations were carried out at the CCSD(T)-F12b level of theory with a cc-pVQZ-PP-F12 basis set for the arsenic atom, and cc-pVQZ-F12 for the hydrogens. Due to the sevenfold increase in computational demand of the DMS over the PES, dipole moments were calculated on a reduced grid of 10,000 points, generated by randomly sampling our PES grid.

5.2 Analytic representation

In a similar spirit to our PES, the ab initio DMS was expressed analytically using the symmetrized molecular bond (SMB) representation31. In this representation the electronically averaged dipole moment $\bar{\bf{\mu}}$ is constructed as symmetrized projections onto the molecular bonds with the dipole moment components ($\overline{\mu}_{A_{1}}$,$\overline{\mu}_{E_{a}}$,$\overline{\mu}_{E_{b}}$) in the molecule fixed axis system given by 4th order polynomial expansions

[TABLE]

where $\Gamma=A_{1},E_{a}$ and $E_{b}$ are the irreducible components of C3v,

[TABLE]

[TABLE]

$\mu^{\Gamma(s)}_{ij...}$ are the expansion parameters, $\Delta r_{k}=r_{k}-r_{\rm eq}$ and $\rho$ is the same as in Eq. (1). The dipole moment components ($\overline{\mu}_{E_{a}},\overline{\mu}_{E_{b}}$) are transformed as linear combinations of each other by the C3v group operations and so transform together as $E$-symmetry. For this reason the parameters $(\mu^{E_{A}(s)}_{ij...},\mu^{E_{b}(s)}_{ij...})$ must be fit together and $\mu^{A_{1}(s)}_{ij...}$ are fitted separately. For an extensive discussion on the SMB representation of the dipole moment function the reader is directed to [ 31].

The final fit required 261 parameters and reproduced the ab initio data with an RMS difference of 0.0008 Debye for energies up to 12,000 cm*-1*, which is comparable to the level of numerical noise in the finite differences procedure. The DMS expansion parameter set and a Fortran 90 routine to construct the DMS is included in the supplementary material.

6 Results

6.1 Structural parameters

Table 3 shows the various structural parameters of AsH3 computed at different levels of theory, compared to those of our refined PES and those derived from experiment. Ab initio calculations of the equilibrium values of $r$ and $\alpha$ were performed using the geometry optimisation procedure in MOLPRO. Both VQZ-PP-F12 and AVQZ-DK level calculations are seen to somewhat overestimate $r_{\rm eq}$ and $\alpha_{\rm eq}$ when compared experiment, a feature that is exacerbated by the exclusion of relativistic effects altogether ($r_{\rm eq}^{\rm AVQZ}=1.52375$ Å, $\alpha_{\rm eq}^{\rm AVQZ}=92.5553^{\circ}$ and $r_{\rm eq}^{\rm AV5Z}=1.523653$ Å, $\alpha_{\rm eq}^{\rm AV5Z}=92.54910^{\circ}$). As expected, the effect of our equilibrium geometry adjustment results in equilibrium bond lengths and angles much closer to that of experiment 65. This is reflected in the good agreement between our purely rotational energies and spin-statistics averaged hyperfine resolved rotational energies of Tarrago et al. 62 (see Table 4). There are small systematic residuals as large as 0.01 cm*-1*, suggesting our treatment of the rotational motion could be improved by further tweaking the equilibrium parameters. However, doing so would undoubtedly spoil the vibrational accuracy so we decided against it.

As yet, the inversion barrier height $\Delta E$(barrier) of AsH3 remains unmeasured. The previous highest-level predictions are those by Schwerdtfeger et al. 66 in 1992, who calculated a value of 13079.3 cm*-1* at the Moeller-Plesset (MP2) level of theory. This is somewhat lower than our CCSD(T) values of just over 14 000 cm*-1*, shown in Table 3. The minimum energy path over the barrier reduces the As-H bond lengths to their so-called saddle-point value $r_{\mathrm{SP}}$ at planar geometry. Of this, the predicted value of 1.457 Å  by Schwerdtfeger et al. is in reasonable agreement with our own (see Table 3). For comparison, the NH3 barrier height is measured to be 1786.8 cm*-1* occurring for $r_{\mathrm{SP}}=0.99460$ Å 67, and for PH3 the calculated values of Sousa-Silva et al. 39 are currently the most reliable, predicting a value of 11 130 cm*-1* at 1.3611 Å. Whereas the NH3 inversion splitting is well known to be $\approx 0.79$ cm*-1* for the ground vibrational state, it has been predicted but not observed in PH3 39, 68 and so it is unlikely to be observed in AsH3 for some time.

6.2 Rovibrational energies

Rovibrational energy level calculations were performed up to $J=30$ using the AsH3-CYT18 PES in conjunction with the nuclear motion program TROVE. Model input parameters were kept the same as reported in section 3, including our P${}_{\mathrm{max}}=14$ vibrational basis. With a basis set of this size the vibrational Hamiltonian E-symmetry block has 2571 roots. Therefore, given the 2$J$+1 multiplication factor for rotationally exited states, it was necessary to perform additional basis set truncations to reduce computational cost. Firstly, our purely vibrational energies $E^{i}_{\rm vib}$ were truncated at 12 000 cm*-1*. These, upon multiplication with rigid symmetric rotor wavefunctions, form the basis for our full rovibrational calculation, which we term the $(J=0)$-contracted basis. Our second truncation, performed only once $J$ exceeded 21, is therefore to remove all $(J=0)$-contracted eigenfunctions with energy greater than $E^{i}_{\rm vib}+E^{i}_{\rm rotor}=$16 000 cm*-1*, where $E^{i}_{\rm rotor}$ are eigenvalues of a symmetric rigid rotor.

Our complete list of calculated energies is available from the ExoMol website (www.exomol.com), along with associated local mode quantum labels $(n_{1},n_{2},n_{3},n_{4},n_{5},n_{6},\Gamma_{\mathrm{vib}},J,K,\Gamma_{\mathrm{rot}},\Gamma_{\mathrm{tot}})$. Here $(n_{1},n_{2},n_{3})$ are stretching quantum numbers, $(n_{4},n_{5},n_{6})$ are bending quantum numbers, $K$ is the projection of the total rotational angular momentum $J$ onto the molecular axis of symmetry, and $(\Gamma_{\mathrm{vib}},\Gamma_{\mathrm{rot}},\Gamma_{\mathrm{tot}})$ are the vibrational, rotational and total symmetry in C3v. The local mode vibrational quantum numbers can be converted to the normal mode representation using symmetry rules (see Down et al. 69), under the assumption that the total number of stretching and bending quanta are conserved between representations. Included in the supplementary material is a list of calculated vibrational states that have been converted to the normal mode representation for $n_{1}+n_{2}+n_{3}\leq 4$ and $n_{4}+n_{5}+n_{6}\leq 4$. This covers all strong bands under 7000 cm*-1*, and should aid any future labelling of AsH3 spectra.

Table 5 compares the calculated $J=0$ term values under 7000 cm*-1*, computed using our refined PES, AsH3-CYT18, to the experimentally observed values 16, 58. Vibrational labels above 5000 cm*-1* are tentative. The refined PES reproduces empirical energies with a marked improvement over the ab initio surface (see Table 2). All bands included in the refinement, except for the $2\nu_{4}^{l=2}$ band, display sub-wavenumber accuracy. Based on our energy residuals for $J=1-6$ states belonging to the $2\nu_{4}^{l=2}$ band, which fall within $\pm$1.0 cm*-1* of experiment, we strongly suspect the empirical band origin20 of 2003.483 cm*-1* is incorrect. From a comparison of our $J=1,2$ energy residuals, we expect the true value to be closer to 1997.5 cm*-1*. Interestingly, for the bands at 3000 and 5000 cm*-1* measured by Halonen et al. 58, all calculated $J=0$ term values, except for the $2\nu_{1}+\nu_{2}$ and $\nu_{1}+\nu_{2}+\nu_{3}$ bands, fall within a few wavenumbers of experiment despite being omitted from the refinement. This illustrates the interpolative power of the refinement, and suggests that even bands not yet observed experimentally should be predicted with reasonable accuracy by our refined PES. Alternative matches for the 5057 cm*-1* bands within our energies list would be the $\nu_{1}+3\nu_{4}^{1}$ (predicted at 5052.561 cm*-1*) and $\nu_{3}+3\nu_{4}^{1}$ ($A_{2}$) (5052.758 cm*-1*) bands. However, considering that these bands are not predicted to be observable at room temperature, the discrepancies are more likely due to resonance interactions that are not well modelled in our PES.

The residual differences between our calculated $J=0$ term values and those of experiment can be removed from the final line list by utilising an empirical basis set correction (EBSC) 31, whereby our calculated band centres are simply replaced by the corresponding experimental values. We employed the EBSC for all experimentally known bands taken from 16, except the suspicious $2\nu_{4}^{l=2}$ band. Figure 2 displays the difference between our calculated energies and those derived from experiment for states with $J\leq 6$ taken from Refs. 62, 18, 19, 20, 21, 22 after employing the EBSC. The corresponding root-mean-square errors ($\sigma_{\rm rms}^{\rm ebsc}$), split by vibrational band, are shown in table 5. Although there is some deterioration in quality with $J$, this is slow and systematic in most cases, reassuring us that our calculations can safely be extended to higher rotational excitations. Agreement for the $2\nu_{2}$ and stretching bands is particularly pleasing, and all calculated $J=1-6$ energies, bar those belonging to the $2\nu_{3}^{2}$ band, are calculated to within $\pm 0.2$ cm*-1* of the experimental values. Judging by the systematic offset of the $2\nu_{3}^{2}$ band in figure 2, the experimental band center used in the EBSC is most likely $\approx 0.2$ cm*-1* lower than the true value. Slightly larger $\sigma_{\rm rms}^{\rm ebsc}$ values are observed for the $2\nu_{4}^{0}$, $2\nu_{4}^{2}$ and $\nu_{2}+\nu_{4}$ bands. Whereas the $2\nu_{4}^{2}$ and $\nu_{2}+\nu_{4}$ bands display clear $J-K$ dependencies, it is difficult to discern any such trends for the $2\nu_{4}^{0}$ band, which was omitted from the refinement altogether. Two possible reasons for this are either corrupt experimental data, or perturbation interactions due to nearby states that are not correctly represented by our PES. Even so, the 0.207 cm*-1* root-mean-square error is very reasonable.

6.3 Line intensity predictions

To simulate absolute absorption intensities we use the expression

[TABLE]

where $\Phi_{int}^{\prime}$ and $\Phi_{int}^{\prime\prime}$ are the upper and lower state wavefunctions respectively, that correspond to energies $E^{\prime}$ and $E^{\prime\prime}$. $\nu_{if}$ is the transition frequency, $\mu_{A}$ is the electric dipole moment along the $A=X,Y,Z$ axis, T is the absolute temperature and $Q(T)$ the partition function given by

[TABLE]

In the above equation, $E_{w}$ is the energy and $g_{w}$ is the total degeneracy of state $w$. Note that we are using the same symbol $E$ for the energies (Joule) and term values or wavenumbers (cm*-1*), whereas conventionally the latter should be $\tilde{E}$.

The nuclear spin statistical weights for AsH3 are {16,16,16} for states of {$A_{1}$,$A_{2}$,$E$} symmetry, and so $g_{w}=16(2J_{w}+1)$. No calculated or experimentally derived values of the partition function could be found in the literature, so we provide partition function values in the supplementary material, for temperatures ranging from 10 to 500 K in intervals of 10 K. Fig. 3 illustrates the convergence of $Q$ as the rotational basis is increased from including only $J=0$ states ($J_{\rm max}=0$), to all computed states with $J\leq 30$ ($J_{\rm max}=30$). In reality there will be additional contributions from our vibrational basis (P${}_{\rm max}=14$) and PES, although these are difficult to quantify. The room temperature partition function was calculated to be $Q(T=296)=8250.2801$ using $J_{\rm max}=30$, which we estimate to be better than 99% converged.

Line list calculations were performed using the AsH3-CYT18 PES and the cc-pVQZ-PP-F12 DMS detailed in section 5. Transitions involve states with energies up to 10 500 cm*-1*, rotational excitation up to $J=30$, and a maximum lower state energy of 3500 cm*-1*. The final line list consists of 3.6 million absorption lines in the range $0-7000$ cm*-1* with intensity greater than 1$\times 10^{-28}$ cm*-1*/(molecule cm*-2*) at 296 K. An overview is presented in figure 4. It is available to download from the ExoMol website (www.exomol.com), where it is provided in the ExoMol format 70. Summarily, this consists of a .states file which contains the complete list of rovibrational states with associated energies and quantum numbers, and a .trans file which contains the complete list of transitions, each identified by an upper and lower state index (in the .states file), Einstein A-coefficient, and transition wavenumber.

Several sources of experimental absorption data exist for AsH3. In the following paragraphs our intensity calculations are validated by comparison with only the most recent and reliable sources. For the first test of our absolute intensities we compare our calculated band intensities with those obtained by Zheng et al. 71, shown in Table 6. Zheng et al. produced a three-dimensional DMS based on of density functional theory calculations, and compared the resulting absolute vibrational band intensities to the values obtained by direct integration of absorbance spectra, which they provide with 20$-$40% estimated uncertainty. Due to multiple bands overlapping only the combined intensity of bands with the same local mode quanta are presented in some cases. For the $\nu_{1}$ and $\nu_{3}$ fundamentals we compare well with experiment, reproducing the observed values within 2% and 14% respectively. Zheng et al. only provide the measured intensity of the sum of the $2\nu_{1}$ and $\nu_{1}+\nu_{3}$ bands, for which we are stronger by 17%. No measurements of the weaker $2\nu_{3}^{l=0}$ and $2\nu_{3}^{l=2}$ bands are given, most likely due to difficulties resolving them without accurate theoretical line positions. Finally, for the three-quanta stretches, our calculated intensities are typically 2-3 times weaker than the measured values. However, it is difficult to estimate the reliability of these measurements, given the recorded spectrum is only medium-resolution ($\Delta\nu=0.2$ cm*-1*) and the bands are weak.

Dana et al. 15 measured absolute intensities of 387 lines belonging to the $\nu_{1}$ and $\nu_{3}$ bands. Their line measurements range from 2010$-$2235 cm*-1* although they make no attempt to measure the $Q-$branch from 2110$-$2140 cm*-1*, presumably owing to the density of lines. Table 7 compares our calculated line positions and intensities to the experimentally measured values for 14 randomly selected strong lines measured by Dana et al. In all cases our calculated intensity values are within $\pm$10% of experiment, although there is a slight tendency to be higher. Nevertheless, this is reassuring given our 14% discrepancy with the $\nu_{3}$ band as measured by Zheng et al. 71.

The PNNL database provides a composite spectrum of pure AsH3 up to 6500 cm*-1* measured at 5, 25 and 50*∘C. For comparison, we generated synthetic T=298.15 K spectra using a $J=0-30$ line list convoluted with a Voigt profile with half-width at half-maxima (HWHM) of 0.09 cm-1*. Although linewidths are well known to depend upon the upper and lower state quantum numbers, the strongest dependency being $J$ and $K$, as far as we know no such data exists for AsH3, and we found the value 0.09 cm*-1* reasonably approximated the PNNL linewidths on average. To convert the PNNL absorbance spectra to cm2/molecule a multiplication factor of 9.28697$\times 10^{-16}$ is necessary.

Figure 5 shows an overview of our synthetic spectrum compared to PNNL, with the key absorption features expanded in Figs. 6-15. Qualitative agreement is very good, particularly for the $\nu_{1}$/$\nu_{3}$ (see Figs. 8 and 9) and $\nu_{2}$/$\nu_{4}$ (see Figs. 6 and 7) fundamentals, and the $2\nu_{1}$/$\nu_{1}+\nu_{3}$/$2\nu_{3}^{0}$/$2\nu_{3}^{2}$ band system (shown Figs. 12 and 13). Note that despite only including the $\nu_{2}$ and $\nu_{4}$ band centres in the refinement, their rotational structures are reproduced well. In the 2920$-$3260 cm*-1* region (shown Figs.10 and 11) the dominant sources of opacity are predicted to be the strong $\nu_{1}+\nu_{4}$ and $\nu_{3}+\nu_{4}$ (A1) (calculated band centre 3119.400 cm*-1*) bands, and the slightly weaker $\nu_{2}+\nu_{3}$ (3023.706 cm*-1*) and $\nu_{1}+\nu_{2}$ bands. Considering that no associated experimental energies were included in the refinement, the level of agreement is satisfying. Above 5000 cm*-1* most absorption features are lost in the PNNL background noise; only the $2\nu_{1}+\nu_{4}$, $\nu_{2}+2\nu_{3}^{0}$ and $\nu_{1}+\nu_{3}+\nu_{4}$ bands (our labelling) are clearly visible between 5000–5250 cm*-1* (see Fig. 14). There is a tenuous absorption bump in PNNL at 5050 cm*-1* for which we appear to be offset by roughly 15 cm*-1*, confirming our discrepancies with the $2\nu_{1}+\nu_{2}$ and $\nu_{1}+\nu_{2}+\nu_{3}$ band centres measured by Halonen et al. 58 (see table 5). In the $6000-6400$ cm*-1* region (see Fig. 15) the salient feature is the $3\nu_{1}$ and $2\nu_{1}+\nu_{3}$ Q-branch at 6135 cm*-1*, for which we clearly underestimate the intensity. Although in line with our comparisons with Zheng et al. 71 (see table 6) it is difficult to quantify this, or indeed draw any conclusions regarding the weaker $\nu_{1}+2\nu_{3}$/$3\nu_{3}$ bands, without additional high-resolution experimental data.

The largest source of error in our intensity calculations will undoubtedly be the DMS. To improve on the CCSD(T)-F12b/cc-pVQZ-PP-F12 method by Hill et al. 46, large CCSD(T)/aug-cc-pwCVnZ-DK ($n=4,5$) calculations would likely be necessary (for an example, including additional post-CCSD(T) corrections see Ref.  72), which are currently computationally unmanageable for a full 6D surface. Secondarily, the PES quality must be considered. Line intensities are inexorably connected to the PES through the wavefunctions in Eq. (20), and accurate modelling of intensity transfer between lines (so-called ‘intensity stealing’) relies on the correct representation of rotation-vibration resonances in the PES. Therefore, from a nuclear motion point-of-view, high-resolution measurements complete with line intensities and quantum assignments, particularly for the $800-1200$, $2900-3300$ and $5000-5300$ cm*-1* regions, would be most beneficial for future modelling.

7 Conclusion

We have produced the first full-dimensional PES and DMS for the arsine molecule. Both PES and DMS were computed at the CCSD(T)-F12b/cc-pVQZ-PP-F12 level of theory, with implicit treatment of scalar relativistic effects via replacement of 10 core electrons with a relativistic pseudopotential. A comparison with standard CCSD(T)/aug-cc-pVQZ-DK based calculations employing the DKH8 Hamiltonian, showed that CCSD(T)-F12b/cc-pVQZ-PP-F12 level theory resulted in significantly more accurate nuclear motion calculations.

Geometry optimisation and empirical adjustment of harmonic and certain cubic terms in the pVQZ-PP-F12 PES resulted in $J=1-6$ rotational energies with a root-mean-square error of 0.0055 cm*-1*, and vibrational term values accurate to within 1 cm*-1* for all reliably known experimental band centres under 6400 cm*-1*. Utilising the empirical basis set correction scheme, 578 experimentally derived ($J=1-6$) rovibrational energies are reproduced with an RMS of 0.122 cm*-1*. Vibrational term value comparisons with eight approximately known band centres showed that six agreed within 3.5 cm*-1* despite being omitted from the refinement. The remaining two displayed $\sim 16$ cm*-1* discrepancies, most likely due to overlooked resonances.

Rotational-vibrational line intensity calculations were performed using the refined PES and ab initio DMS, in conjunction with variational nuclear motion calculations. The resulting line list, with full quantum assignments, extends to 7000 cm*-1* and is complete up to 300 K. Comparisons with multiple experimental sources show our intensity predictions to be reliable, in particular, good overall agreement with the main absorption features present in PNNL is noted. Our complete line list with quantum assignments is available from the ExoMol website (www.exomol.com) in ExoMol format 70.

As far as we know, arsenic is the heaviest element for which there exists an associated variationally-computed infrared molecular line list. Considering that the quantum chemistry methods employed here are available for most p-block main group elements 46, the outlook for studying similar systems in future is positive.

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

PAC thanks EPSRC for a CASE studentship under grant EP/L504889/1 and Servomex for industrial sponsorship. The authors are grateful to Grant J. Hill for his help with electronic structure calculations. SY and JT thank STFC for support under grant ST/M001334/1. The authors acknowledge the use of the UCL Legion High Performance Computing Facility (Legion@UCL), and associated support services, in the completion of this work, along with 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..