Using block diagonalization to determine dissociating autoionizing states: Application to N2H, and the outlook for SH

We describe our implementation of the block diagonalization method for calculating the potential surfaces necessary to treat dissociative recombination (DR) of electrons with N2H+. Using the methodology we have developed over the past few years, we performed multi-reference, configuration interaction calculations for N2H+ and N2H with a large active space using the GAMESS electronic structure code. We treated both linear and bent geometries of the molecules, with N2 fixed at its equilibrium separation. Because of the strong Rydberg-valence coupling in N2H, it is essential to isolate the appropriate dissociating, autoionizing states. Our procedure requires only modest additional effort beyond the standard methodology. The results indicate that the crossing between the dissociating neutral curve and the initial ion potential is not favorably located for DR, even if the molecule bends. The present calculations thereby confirm our earlier results for linear N2H and reinforce the conclusion that the direct mechanism for DR is likely to be inefficient. We also describe interesting features of our preliminary calculations on SH.


Introduction
In recent work [1][2][3] we have shown that the powerful techniques of quantum chemistry can be adapted to the calculation of dissociating, autoionizing states by implementing the block diagonalization method.This manuscript describes calculations that address the dissociative recombination (DR) of N 2 H + , a process thought to be quite important in the interstellar medium: EPJ Web of Conferences Because nitrogen in the interstellar medium cannot be directly observed, its abundance is inferred indirectly by detecting N 2 H + , which can be produced by the reaction It was long asssumed that the N 2 consumed in the reaction above was almost completely recovered by dissociative recombination (DR).The measurements of Adams et al. [4] in 1990 indicated that the branching ratios in the dissociative recombination process strongly favor the N 2 + H channel.However, in 2004 Geppert et al. [5] reported measurements that the branching ratio to the NH + N channel was dominant.This unexpected result cast doubt on the assumptions mentioned above and motivated additional experimental and theoretical work.
In 2007, Molek et al. [6] reported a new experimental measurement of the branching ratio and concluded that the upper limit for the branching ratio to the NH + N channel was 5%.Also in 2007, Talbi [7] reported a theoretical investigation of the potential surfaces for linear N 2 H and N 2 H + .She concluded that the likely outcome of dissociative recombination would be the N 2 + H channel, with N 2 in the first electronically excited state.There is now general agreement that the branching ratio to the NH + N channel is less than about 5% [8].
Talbi's work [7,9] highlighted the need for further calculations to determine appropriate diabatic potential surfaces for both linear and bent geometries of N 2 H.Our work in 2012 [3] reported detailed calculations for the linear molecule that suggested that the N 2 + H channel was the more important.This manuscript investigates this issue further by considering bent geometries of N 2 H.

Determination of autoionizing states using block diagonalization
A detailed calculation of dissociative recombination must provide potential surfaces for both the initial molecular ion and the dissociating autoionizing state.Calculating the potential surface for the initial molecular ion (N 2 H + in the present case) is straightforward and involves standard techniques of quantum chemistry.However, calculating the neutral dissociating state is much more difficult because that state is embedded in the continuum of scattering states that correspond to e − + N 2 H + .
Recent work by ourselves [1][2][3]10] and others [11,12] has shown that the block diagonalization method provides a very powerful technique to determine the necessary dissociating, autoionizing states and coupling terms for N 2 H.The principal difficulty that one encounters when trying to calculate the dissociating surface using the standard techniques of electronic structure is that the states of interest are (1) highly excited and (2) strongly mixed with other states.Multiple curve crossings occur, and it can be very difficult to isolate the desired potential surface.
The block diagonalization method [11] provides an effective technique for transforming the results of a standard electronic structure calculation into diabatic potential curves, thereby unraveling the complicated pattern of interactions and identifying the dissociating state.An advantage of the method is that one can perform a conventional electronic structure calculation of the desired size and accuracy and then obtain diabatic curves with comparable reliability.A second advantage is that the numerical effort for the diabatization is quite modest compared to the original calculation (which usually involves diagonalizing a very large, sparse matrix.) The method requires extra effort in the determination of molecular orbitals (mo's).Briefly, one must insure that the variation of the mo's { i } with molecular geometry is small.This result can be achieved by setting up the calculation in such a way that the adiabatic energies are invariant under a rotation of certain molecular orbitals, and then selecting an appropriate rotation at each geometry to transform the orbitals into a standard form.The details of the methodology are available in the literature [1][2][3]11].The final result is that large Hamiltonian matrix of size N × N determined by the CI calculation is brought into block diagonal form, as shown in Fig. 1.The small block on the upper left of the final matrix is the desired diabatic matrix of size N × N , where N N .An important numerical consideration is that after the desired energies and eigenvectors of the large CI matrix have been found by standard electronic structure techniques, the additional calculations require only operations on matrices of order N .The diagonal elements of the diabatic matrix correspond to energies for the valence and Rydberg states.One can easily identify the dissociating states.The off-diagonal elements provide Rydbergvalence coupling terms, and one can scale these to estimate the coupling V el of the autoionizing valence states to the background continuum.The result is V el ∼ = (n * ) 1.5  Rydberg |H elec | dissoc , where n * is the effective quantum number of the Rydberg state determined by its binding energy relative to the parent ion.The autoioinizing width is 2 |V el | 2 , where is the density of states; 2 = 4/k, and k is the wave number of the continuum electron.

Electronic structure calculations for N 2 H
The electronic structure calculations reported here were performed using the GAMESS code [13].The basis set included a 6 − 311G(d, p) on each N plus diffuse functions and a 10s4p1d basis on H; there were 69 total basis functions.Using the notation of Talbi [7] for the molecular orbitals, the orbital occupancy of the ground state of N 2 H + is where the vertical bar (|) signifies the highest orbital occupied in the ion core.The orbital occupancy of the lowest N 2 H dissociating state is The Rydberg states of N 2 H have one electron in a highly excited orbital outside of the ion core: we use a light gray to indicate that orbitals such as (NH * ) are not occupied.The excited dissociating valence states of N 2 H have an unoccupied core orbital: x * 1 y .The dissociating state important for DR is a linear combination of states corresponding to the two orbital occupancies given above.
An important requirement for the block diagonalization method is that the molecular orbitals vary smoothly with geometry.For most of the calculations reported here, it was sufficient to use natural orbitals determined by a preliminary CI calculation on the lowest N 2 H dissociating state.Energy (E H ) Figure 2 illustrates the smooth variation of the orbitals we used.These orbitals can be determined by the GAMESS code.In a few cases, such as small values of R d , it was necessary to define reference orbitals explicitly as described in [3].

03003-p.3
Figure 3 shows the standard adiabatic curves obtained by GAMESS calculations performed in the linear geometry and also the corresponding diabatic curves.The adiabatic curves exhibit multiple curve crossings; diabatization allows us to identify the appropriate dissociation channels.These calculations were reported earlier [3].The diabatic curves reveal that the crossing between the second repulsive state of N 2 H and the ion is slightly too high to allow for an efficient direct mechanism for the DR of cold N 2 H + .This result led us to investigate whether the crossing might be more favorable for bent geometries of the ion.

DR2013
For the bent molecule, we report new multi-reference, configuration interaction (MRCI) calculations using the orbitals discussed above.The Calculations are somewhat different because of the loss of linear symmetry.We tested two different active spaces for our calculations of the bent molecule.The calculations with the Small Active Space (SAS) have the following structure: Virtual Space : (NN * )(Ryd1)(Ryd2)(Ryd3) . . .
We included single and double excitations from the active space, leading to 449,892 configurations.The calculations were performed using the C s point group, and the potential surfaces reported are of A symmetry.We typically calculated 25-85 eigenvalues and eigenvectors of the electronic Hamiltonian, and the jobs took on the order of 5-24 hours to complete on XSEDE supercomputers [14].
We also performed a smaller number of calculations with a larger active space, denoted the Big Active Space (BAS): Virtual Space : (NN * ) . . .
The calculations with the BAS include Rydberg orbitals in the active space and, as expected, produce more accurate Rydberg potentials.However, the MRCI calculations with the BAS (including single and double excitations) using the C s point group have 3,922,296 CSFs (Configuration State Functions).These jobs are very large and many roots are required (over 70 for values of R d ≤ 1.750 Å).Many of these jobs take between five and 12 days to complete running on XSEDE supercomputers.Fortunately, comparison of the SAS and BAS calculations indicates that the smaller calculations are adequate for calculating the necessary dissociating curves.The comparison is presented in Fig. 4. Note that both calculations shown are s.We primarily investigated the case in which DR leads to breaking the NH bond, so we froze the NN bond length at the equilibrium value appropriate for N 2 H (taken to be R NN = 1.12 Å).For this case it is convenient to use Jacobi corrdinates, which are illustrated in Fig. 4. R d is the distance from the H atom to the center of mass of the N 2 , and is the bending angle ( = 0 corresponds to a linear molecule).One can then expand the potential surface using Legendre polynomials of even order: We typically performed calculations for 8-11 different angles at each value of R d .The calculated points are shown in Fig. 4. We have plotted the points as functions of sin 2 = 1 − cos 2 in order to facilitate least squares fitting.A Legendre polynomial expansion with only even terms for = 0, . . ., 2N is an N th degree polynomial in cos 2 or, equivalently, in sin 2 .By plotting the calculated points vs. sin 2 we could assess the quality of fit using standard polynomial fitting routines.The curves shown in Fig. 4 are typically five-term polynomials (in sin 2 ).It is clear from these curves that an exact fit to the calculated points would have additional undulations.We prefer to have a somewhat smoother fit, because the additional variation in the calculated points is likely due to the uncertainty in selecting the optimum N and the most important configurations for each geometry.Because of the many curve crossings, the set of configurations that make the most important contribution     is constantly changing.Note that the points most distant from the fitted curve tend to be at small values of R d , where the potential is the steepest and small uncertainties are unlikely to be important.
An alternate way to view the diabatic potentials is shown in Fig. 5.One can see that the crossing between the ion and the dissociating curve does not become more favorable as the molecule bends.

Preliminary calculations for the DR of SH +
SH + is observed in various objects of the interstellar medium [15][16][17].Although little is known concerning the interstellar chemistry of this ion, its abundance suggests that the destruction processes for SH + are not efficient.One of the reactions that could destroy it is the electronic dissociative recombination e − + SH + → S + H. (6) 03003-p.6

DR2013
The purpose of the present work is to study this process with the same methodology we applied to the DR of N 2 H + .The ground electronic state of SH + is a 3 − state ( 3 A 2 using the C 2v point group).Its electronic configuration can be written as follows: 1s 2 2s 2 2p 6 3s 2 (SH) 2 1   x 1 y .When an electron recombines with SH + , an excited doublet or quartet state of SH can be formed.We discuss below the SH doublet states.
The doublet states we have considered include the SH Rydberg state with the electronic configuration 1s 2 2s 2 2p 6 3s 2 (SH) 2 1   x 1 y 4s 1 , which is of − (A 2 ) symmetry and correlates with the Rydberg atomic state S(3s 2 3p 3 4s 3 S o ) + H(1s).Other doublet states are the Rydberg states of electronic configuration 1s 2 2s 2 2p 6 3s 2 (SH) 2 1   x 1 y 4p 1 , which can be of − (A 2 ) or (B 1 or B 2 ) symmetry and correlate to the Rydberg state S(3s 2 3p 3 4p 3 P ) + H(1s).This n = 4 SH state is the lowest state in the Rydberg series that converges to the lowest 3 − ionic state of SH + .There is another series of Rydberg states of higher energy, which also converges to the lowest ionic SH + state.This series corresponds to the occupation of the 3d orbitals of SH + , that is, to states of electronic configuration 1s 2 2s 2 2p 6 3s 2 (SH) 2 1   x 1 y 3d, which can be of − (A 2 ), (B 1 or B 2 ) or (A 1 or A 2 ) symmetry.These SH doublet states correlate to the d atomic state S(3s 2 3p 3 3d 3 D o ) + H (1s).
The recombination of SH + with an electron can also lead to the valence states of SH.The lowest one is a state (B 1 or B 2 ) with electronic configuration 1s 2 2s 2 2p 6 3s 2 (SH) 2 3   that correlates with the lowest S(3p 4  where (SH SH * ) 5 means that five electrons are distributed into the three specified orbitals.These configurations lead either to the next higher excited valence states of SH [the 2 + (A 1 ), 2 (A 1 or A 2 ) and 2 (B 1 or B 2 )] and correlate with the S(3p 4 1 D) + H(1s) limit or to the 2 + (A 1 ) state and correlate with the S(3p 4 1 S) + H(1s) limit.
To calculate the potential energy surface governing the dissociation of the SH + ion in its ground electronic state, one must consider at least the ground state of SH, its valence excited states listed above, and the n = 4 and the 3d SH Rydberg states also listed above.In other words, one must use the four lowest − states, the three lowest states, the two lowest states, and the two lowest + states of SH.As a starting point we have set up our ab initio calculations to include the lowest n = 4 Rydberg states.For that, the active space of our MRCI calculation was set up in order to accommodate the 3s, (SH), x , y , (SH * ), 4s, 4p z , 4p y , and the 4p x orbitals of SH.The 1s, 2s, and 2p orbitals of S were frozen.With this active space, using the C 2v point group, and including all single and double excitations, our MRCI calculations involved 3,114,864-3,130,904 CSF's (Configuration State Functions), depending on the state symmetry.For accurate energy calculations we used the Dunning-type correlation consistent basis sets augmented by an n = 4 diffuse function on S and an n = 2 diffuse function on H, that is, the aug-cc-PvTZ basis set as implemented in the GAMESS code used the present calculations.This basis set ensured n = 4 molecular Rydberg orbitals of localized character all along the dissociating coordinate, allowing for their easy identification all along the dissociation path.
Figure 6 shows our most recent calculations of the lowest SH doublet states.Excited SH states of the same symmetry interact strongly, particularly the lowest two excited 2 states.The analysis of the corresponding CI wave functions reveals strong mixing between Rydberg and valence states.
The goal of this preliminary work was to validate the design of our ab initio calculations, and we have obtained very satisfactory agreement of our asymptotic limits with the NIST atomic spectra database of S, as well as with previous calculations of the excited states of SH [18].Our work has also revealed that a thorough treatment of the DR processes of ground state SH + and a diabatization of the states concerned will require the consideration of additional higher Rydberg states, such as the 3d series and also the Rydberg series converging to the first excited 1 state of the ion.These calculations will require extensive ab initio calculations; work along these lines is in progress.

Concluding remarks
We have described large scale electronic structure calculations for the ground state of N 2 H + and for several excited valence and Rydberg states of N 2 H, including the dissociating state important for dissociative recombination.For N 2 H, we have shown that bending the N 2 H + ion does not improve the crossing between the second repulsive state of N 2 H and the ion potential.We therefore do not not expect the direct mechanism for DR to be efficient for cold N 2 H + .We are now investigating the indirect mechanism.
We have also reported preliminary results for SH + and SH; the methodology of our ab initio calculations calculations has been validated by comparison with experimental energy levels and previous ab initio calculations.

Figure 1 .
Figure 1.Schematic illustration of block diagonalization.The N × N block is the desired diabatic matrix.The larger block of size N − N need not be explicitly determined.

Figure 2 .
Figure 2. Variation of the molecular orbitals with molecular geometry.Red indicates that the orbital has a positive sign; blue is negative.The top row shows the antibonding NH * orbital for R d = 1.61 Å, and = 0, 30, 60, and 90 degrees.(R d and are the Jacobi coordinates defined in Fig. 4.) The second row shows the same orbtial for R d = 2.46 Å, at the same angles.The third row shows the A bonding orbital for R d = 1.61 Å, and = 0, 30, 60, and 90 degrees.The smooth changes in the molecular orbitals ensure that the block diagonalization process gives reliable diabatic potential curves.

Figure 3 .
Figure 3. Potentials for linear N 2 H.The left panel shows the adiabatic potentials determined from the GAMESS electronic structure calculations without any additional analysis.The right panel shows diabatic potentials determined using the block diagonalization method.A vibrational wave function for the ion is also shown.

Figure 4 .
Figure 4.The left panel illustrates Jacobi coordinates.The right panel shows the angular dependence of dissociating potential with the big active space (BAS) and with the small active space (SAS).These ab initio diabatic results have not been adjusted in any way.Calculated points and a five-term polynomial fit are shown.

Figure 5 .
Figure 5. Diabatic potential curves for fixed and R NN .

Figure 6 .
Figure 6.Preliminary calculations for SH + and doublet states of SH.
3P ) + H(1s) limit.Correlating with this same limit is the lowest