D matrix analysis of the Renner-Teller effect: An accurate three-state diabatization for NH2

Some time ago we published our first article on the Renner-Teller (RT) model to treat the electronic interaction for a triatomic molecule [J. Chem. Phys. 124, 081106 (2006)]. The main purpose of that Communication was to suggest considering the RT phenomenon as a topological effect, just like the Jahn-Teller phenomenon. However, whereas in the first publication we just summarized a few basic features to support that idea, here in the present article, we extend the topological approach and show that all the expected features that characterize a three (multi) state RT-type'3 system of a triatomic molecule can be studied and analyzed within the framework of that approach. This, among other things, enables us to employ the topological D matrix [Phys. Rev. A 62, 032506 (2000)] to determine, a priori, under what conditions a three-state system can be diabatized. The theoretical presentation is accompanied by a detailed numerical study as carried out for the HNH system. The D-matrix analysis shows that the two original electronic states 2A1 and 2B1 (evolving from the collinear degenerate Pi doublet), frequently used to study this Renner-Teller-type system, are insufficient for diabatization. This is true, in particular, for the stable ground-state configurations of the HNH molecule. However, by including just one additional electronic state--a B state (originating from a collinear Sigma state)--it is found that a rigorous, meaningful three-state diabatization can be carried out for large regions of configuration space, particularly for those, near the stable configuration of NH2. This opens the way for an accurate study of this important molecule even where the electronic angular momentum deviates significantly from an integer value.


I. INTRODUCTION
In 1934 Renner published a detailed study of a linear polyatomic molecule, which possesses a single adiabatic potential energy surface but then, upon bending, splits into two surfaces. 1 As long as a ͑triatomic͒ molecule is limited to its linear configuration we expect the electronic orbital angular momentum quantum number ⌳ប ͑defined with respect to the molecular axis͒ to serve as a good quantum number ͑i.e., to be conserved͒. Moving slightly away from linearity introduces another quantum number, namely, the nuclear angular momentum quantum number, ᐉប. As long as the deviation is small enough both ⌳ and ᐉ serve as good quantum numbers. [1][2][3][4][5] This deviation from collinearity has one additional effect on the molecular structure. All electronic states ͑with a given energy͒, except for the ⌺ state, split into two states with different energies. Renner was the first to study this splitting and the corresponding physical implications ͑for a ⌸ state͒ and therefore the whole phenomenon, which is also well supported by numerous experiments, bears his name.
Introducing the radial coordinate q ͑which measures the deviation from collinearity͒, Renner assumed that for q Ͼ 0 the ⌳ state splits into two states, i.e., E ⌳ ± ͑q͒, where ͑q͒ = 0 q 2 . Renner showed that as long as 0 is small enough ⌳ and ᐉ are, as mentioned earlier, good quantum numbers, but once it becomes large the two quantum numbers lose their relevance as such, and K, defined as K = ⌳ ± ᐉ, becomes the good quantum number. [1][2][3][4][5] In this way the spectrum of a triatomic molecule which for a weak coupling is characterized by the quantum numbers ͑v ͑±͒ , ⌳, ᐉ = ±1͒ turns, for a strong coupling, into a new spectrum of states, characterized by the quantum numbers ͑v ͑±͒ , K = ⌳ ±1͒.

͑1͒
Here, ⌳±ᐉ ͑s e ͉ s͒ are the two electronic eigenfunctions, s e and s are the electronic and the nuclear coordinates, respectively, and is the corresponding ͑nuclear͒ polar coordinate associated with ᐉ. This kind of interaction term is also expected from the Born-Oppenheimer ͑BO͒ treatment. 14,15 Comment: It is important to realize that is not a regular internal body-fixed coordinate but an out-of-plane coordinate associated with a rigid rotation of the molecule around an axis located in the plane formed by the three atoms ͓see, i.e., Fig. 1.
Applying his model, Renner assumed that p ⌳ ±ᐉ ͑s͒ is a constant independent of s and established that in the vicinity of the collinear axis p ⌳ ±ᐉ ͑s͒ is equal to ⌳, which is an integer. In case of a ⌸ state we have ⌳ = 1 and therefore In the earlier years until the end of the 1950s this value for p ⌳ ±ᐉ ͑s͒ was employed in numerous applications. Later, efforts were made to calculate this value as a function of q and to apply the actual value ͑which was found to differ from 1 or from any other integer͒. [6][7][8][9][10][11][12][13] However, this ͑noninteger͒ calculated value cannot be used in a straightforward way due to difficulties related to the diabatic potentials, which become multivalued. 13 To avoid this difficulty the value of p ⌳ ±ᐉ ͑s͒ is assumed to be ⌳ = 1 in those expressions that are directly related to the diabatic potentials, whereas in other instances the calculated value is used. 6,7,[9][10][11][12][13] It seems that all the numerical treatments related to the Renner-Teller ͑RT͒ interaction are based on the two-state approximation and in no situation was a three-state diabatization carried out.
Our main interest in the Renner-Teller-type molecules is that due to the vibronic coupling, their electronic state of a given energy splits into two ͑or more͒ states ͑with different energies͒, and that this splitting disappears once q =0, namely, when the molecule becomes linear. This situation is reminiscent of the Jahn-Teller ͑JT͒ interaction 3,16-28 which stems from the fact that two ͑or more͒ electronic states become degenerate. The similarity between the two interactions is even more pronounced because just like the points of degeneracy within the JT framework that arrange themselves along seams, 29͑a͒ the degeneracy points within the RT framework are located along the collinear axis which for all practical purposes can be considered as a seam.
To continue, we switch from Renner's notation for the NACT ͓presented in terms of p ⌳ ±ᐉ ͑s͔͒ to our notation that has been used, now for several years, while treating the JT effect. We consider the NACT, jk ͑s͒ ͑known to be a vectorial matrix element͒, which couples the two adiabatic states j and k, and is given as a function of the coordinates of the nuclei. These NACTs are of the form jk ͑s͒ = ͗ j ͑s e ͉sٌ͉͒ k ͑s e ͉s͒͘, ͑3͒ where j ͑s e ͉ s͒ and k ͑s e ͉ s͒ are the corresponding two ͑elec-tronic͒ eigenfunctions and ٌ is the grad operator.
Recently we published a Communication 30 in which we discussed the RT interaction and how it is related to the corresponding RT NACTs and the degeneracy points. To be more specific we considered the HNH molecule and employed an ab initio package to calculate p ⌳ ±ᐉ ͑s͒ for the case ᐉ = ⌳ = 1. In what follows this particular magnitude is designated 12 ͑s͒ and is written in the form FIG. 1. Ab initio RT nonadiabatic coupling terms, 12 ͑q ͉ z͒ and 23 ͑q ͉ z͒, as calculated for the collinear and planar HNH molecule. The results are presented as a function of q-the vertical distance of the rotating atom from the fixed ͑collinear͒ axis-for different configurations: ͑a͒ the rotating atom is nitrogen and the calculations are done for z = 1.95 a.u. ͑the C 2v symmetry͒; ͑b͒ the rotating atom is nitrogen and the calculations are done for z = 1.0 a.u.; ͑c͒ the rotating atom is hydrogen and the calculations are done for z = 3.9 a.u. ͑-͒ 12 ͑q ͉ z͒; ͑¯͒ 23 ͑q ͉ z͒.

͑4͒
To carry out this calculation we assume that two out of the three atoms are clamped at their positions and only one atom is allowed to rotate around the molecular axis. Consequently 12 ͑s͒ described in terms of three cylindrical coordinates related to this atom, namely, s = ͑ , q , z͒ ͑see Fig. 1͒. Locating the origin at the position of the left hand side clamped hydrogen and assuming the z axis to be along the molecular axis, we get that ͑1͒ the coordinate z measures, for the collinear configuration, the distance of the rotating atom from the origin ͑2͒ the radius q measures its vertical distance from the ͑molecule͒ axis; and ͑3͒ the angle is the corresponding ͑out-of-plane͒ cylindrical angular coordinate around the axis ͓see Fig. 1͔. It is important to emphasize that such calculations were carried out not only for the two ⌸ states ͑designated as states 1 and 2͒ but for any two states included in the numerical study.
Preliminary results along these lines ͑for the two ⌸ states͒ were published in that Communication. 30 For completeness they are presented once again ͑following an improved convergence͒. Thus, we consider three different configurations: two configurations refer to the case that nitrogen is rotating around the fixed HH axis ͓see panels ͑a͒ and ͑b͔͒ and one configuration refers to the case that hydrogen is rotating around the fixed HN axis ͓see panel ͑c͔͒. All calculations are done for the case that at the collinear arrangement the distance R HH , between the two hydrogens, is R HH = 3.9 a.u. In panel ͑a͒ are presented results for which the nitrogen is located along the symmetry line ͑the C 2v configu-ration͒ and therefore z = 1.95 a.u.: in panel ͑b͒ are presented results for the case where the nitrogen is shifted away from the symmetry line, namely, for z = 1.0 a.u. ͑thus, the nitrogen is much closer to one of the hydrogens͒. In panel ͑c͒ we show results for the case where a hydrogen and the ͑single͒ nitrogen are clamped whereas the second hydrogen is allowed to rotate around the axis ͑in this case the HN axis͒.
The two main features to be noted are ͑1͒ the value 12 ͑q ϳ 0.0, z͒ϳ1.0 for all three cases; this implies ͑in "our" language͒ that the RT NACTs, in the slightly bent HNH system, are quantized and that their value is 1 ͑and not, e.g., 1 / 2 like in the cases of JT NACTs͒. This result is expected, as elaborated in the literature. 1-5 ͑2͒ However, it is also noted that 12 ͑q , z͒ decreases monotonically as q increases ͑i.e., as the molecule becomes more and more bent͒; but the rate of decrease is z dependent. This feature too is mentioned in the literature 6,7,[9][10][11][12][13][14] ͓mainly for the C 2v configuration presented in panel ͑a͔͒.
It is well known that the deviation of 12 ͑q , z͒ from 1 leads to multivalued diabatic potentials. [6][7][8][9][10][11][12] However, it is not clear if a rigorous treatment was ever done to resolve the multivaluedness of the diabatic potentials. We intend to show that our theory that supplies the means to rigorously resolve such difficulties in case of the JT interaction, applies also for the RT interaction.
The present article contains the following sections: In Sec. II is given a brief account of the theory as developed for the JT intersections and their NACTs. Certain aspects of this theory are then applied for the RT NACTs. In Sec. III are given numerical results as obtained for the HNH molecule: We present potential energy curves and the RT NACTs as a function of q, for different z values, related to the three lower states of this molecule, namely, the two ⌸ states and the state that evolves from the ⌺ state. A discussion and the conclusions are given in Sec. IV in which we also briefly relate to JT cis located in the same region of configuration space.

A. Background comments
As mentioned earlier the theoretical treatment applies to the NACTs due to both the RT interaction and the JT interaction. The RT interaction, like the JT interaction, can be traced back to the BO treatment 14,15 but what really makes the two types similar is the fact that both are caused by the degeneracy of eigenfunctions at certain points in configuration space. Moreover, in both cases, these points have the tendency to arrange themselves along continuous, usually infinite, lines defined in configuration space. These lines are known, in case of the JT interaction, as seams and we suggest applying the same name to the lines in case of the RT interaction. The main difference between the two kinds is that the JT seams are usually accidental and therefore have to be exposed numerically by trial and error, whereas the RT seams are a result of bending modes and therefore run along well defined geometrical lines such as the axes of triatomic molecules, etc.
The theory we intend to present is based on closed contours that surround the seams. 20͑b͒,29͑b͒ The main outcome of this theory is that it supplies the conditions for a group of states to yield, in a given region in configuration space, single-valued diabatic potentials. This theory is successfully applied in case of the JT interaction and we intend to show that it applies also in case of the RT interaction.
The difference in treating these two kinds of interaction is the positions of the seams with regard to the plane formed by the three atoms ͑A, B, and C͒. In case of the JT interaction the seam intersects that plane at one point ͑or at most, at a few isolated points͒ and therefore the contour that surrounds the seam can be assumed, without detracting from generality, to be fully located in that plane. In case of the RT interaction the seam is always located in the plane and therefore the ͑closed͒ contour that surrounds it has to be located outside that plane. The simplest choice is to assume it to be in a plane perpendicular to the seam ͑see Fig. 1͒, although any other contour ͑that surrounds the seam͒ is expected to yield similar results.
Although the seams are located in different positions ͑with regard to the triatom plane͒ the contours are constructed by clamping two atoms and permitting the third atom, that is atom A, to move freely and form the ͑closed͒ contour. In both cases we assume the contour to be in a plane and therefore if the origin of the system of coordinates is in that plane the position of the free atom is determined by two coordinates ͑usually polar coordinates͒. This description can be easily generalized to more elaborate contours.
In both cases the free atom ͑atom A in our case͒ is used to examine the values of the various NACTs, jk ͑s͒; j , k = ͕1,N͖ where N is the number of states included in the manifold and s is a nuclear coordinate in the region of interest. According to our way of presenting the NACTs the point s is written as s = ͑ , q ͉ z͒ where ͑ , q͒ are the two polar coordinates mentioned earlier and z is an indirect coordinate ͑or describes numerous indirect coordinates͒ to characterize different contours.
In the present study we consider systems with two or three states, thus N =2,3, and treat only more significant NACTs, namely, 12 ͑s͒ and 23 ͑s͒ ͓the third NACT, i.e., 13 ͑s͒, is usually negligibly small and therefore is ignored in our study͔.
Before we start elaborating on the theoretical tools to be applied in this study we mention another issue, namely, the Hilbert subspace. 20͑b͒,29͑b͒ A Hilbert subspace, usually in a given region in configuration space, is defined as a group of states which approximately are not coupled ͑or at most are weakly coupled͒ to states outside the group. Thus if jk ͑s͒ is a NACT that couples states j and k where state j belongs to the subspace and state k does not, then for all practical purposes all its components are negligibly small. The importance of the subspace stems from the fact that the corresponding NACTs ͓see Eq. ͑3͔͒ can form single-valued diabatic potentials. In others words if a group of states does not form a Hilbert subspace the corresponding NACTs will deliver multivalued diabatic potentials which cannot be applied in any further studies.
In the next subsections we briefly discuss the connection between the ability to form single-valued diabatic states and the features of the nonadiabatic coupling matrix ͑NACM͒. This was done on numerous occasions in case of the JT NACM and here, among other things, we emphasize the friendly features of the RT NACM.

B. The NACM and the diabatic potential matrix
Having the BO ͑diagonal͒ adiabatic potential matrix, u͑s͒, the diabatic potential matrix W͑s͒ is obtained following the adiabatic-diabatic transformation ͑ADT͒ matrix A͑s͒: 20͑b͒,29͑c͒,31 The ADT matrix A͑s͒ can be shown to be an orthogonal ͑unitary͒ matrix that fulfills the following first order differential ͑vector͒ equation 31 where ͑s͒ is the NACM with the elements as defined in Eq. where ဧ is the ordering operator, s 0 is the initial point of integration, ⌫ is the contour along which Eq. ͑6͒ is required to be solved, the dot stands for a scalar product, and A͑s 0 , ⌫͒ is the initial value of A͑s͒ on ⌫. In what follows A͑s 0 , ⌫͒ is assumed to be the unit matrix. It is well noticed that the only component of ͑s͒ that affects the above line integral is the tangential component along the contour ⌫.
Another matrix of interest is the topological matrix D͑⌫͒ which is identical to the A matrix but is calculated for a closed contour. 20͑b͒,23͑b͒,23͑c͒,29͑b͒,29͑c͒ Thus It can be shown that in order for the diabatic potential matrix W͑s͒ to be single valued in the region of interest the D matrix, for any chosen closed contour ⌫ in the region, has to be diagonal. Since D ͓just like A͑s͔͒ is unitary its elements are expected to be 20͑b͒,23͑b͒,23͑c͒,29͑b͒,29͑c͒ where j ͑⌫͒; j = ͕1,N͖ are real phases. In case of real eigenfunctions the phases become integer multiples of so that the D-matrix elements are Next we briefly analyze what happens in case ⌫ is chosen to be a circle defined by the position of its center and the relevant radius q. In this situation the ADT matrix can be written as

͑10͒
where we identify ͑1/q͒ ͑ ͉ q͒ as the angular component of . In the same way the topological matrix D takes the form

͑11͒
Earlier we said ͑without proof͒ that the condition for the diabatic potentials to be single valued is the fulfillment of Eqs. ͑9͒. Moreover, it is seen that the single valuedness is solely determined by the NACM. Therefore in order to find out if a group of N states is capable of yielding single-valued diabatic potentials all that is needed is to calculate the corresponding N ϫ N D matrix and see to what extent it is diagonal.
Comment: Equation ͑9b͒ presents two possibilities for the various diagonal D-matrix elements, namely, ±1. These two possibilities are typical for the JT NACTs. However, in case of the RT NACTs we encounter one sign only that can be assumed to be +1.

C. The Renner-Teller topological matrix
Assuming the center of the circular contour ⌫ to be located on the collinear axis, we find that the NACM is independent of the polar coordinate so that the D matrix as given in Eq. ͑10͒ takes the form where the integration over was done trivially because the NACTs, just like the potentials, do not depend on out-ofplane coordinates. To continue we replaced ⌫ by z to specify the circle that is employed ͑this circle is also characterized by q͒. In what follows we discuss two cases, namely, the two-state Hilbert subspace and the three-state Hilbert subspace. Comment: In Eq. ͑12͒ we dropped the subscript because, from now on, we refer only to the angular components of both the NACM and its NACTs and no confusion is expected. Thus and jk stand for and jk , respectively.͔

͑14͒
In order for the D matrix to be diagonal the NACT, 12 ͑q ͉ z͒, has to be 12 ͑q,z͒ = ͭ n; ͑a͒ ͑2n + 1͒/2; ͑b͒ ͮ ͑15͒ where n is an integer. In other words, the two states under consideration form a Hilbert subspace in a region if and only if for each point ͑q , z͒ in that region n is an integer. Later, in the numerical part, we show that case ͑a͒ applies to the RT interaction.

The three-state Hilbert subspace
To treat the three-state case we assume the 3 ϫ 3 matrix to be of the form and we recall that jk ϵ jk ͑q ͉ z͒. It is well noticed that the D ͑3͒ matrix in Eq. ͑17͒ becomes diagonal if and only if = n where n is an integer. It is interesting to mention that the eigenvalues of ͑q ͉ z͒ are ͕i ,−i ,0͖. Since our main concern is the diagonal elements of the D ͑3͒ matrix we calculate them employing the following expressions: Summary: In order for the three states to form single-valued diabatic potentials ͑in the region of interest͒ each of the three matrix elements in Eq. ͑19͒ has to be equal to 1 at each point in that region. Among other things this implies that the eigenvalues of ͑q ͉ z͒ at each such a point are ͕i ,−i ,0͖.

III. NUMERICAL RESULTS
Three kinds of results are discussed in the article: the ͑two͒ tridiagonal elements of the NACM, the corresponding ͑three͒ potential energy surfaces, and the diagonal elements of the topological D matrix. These results are presented as a function of q for the three configurations already given in Fig. 1 ͓details related to Fig. 1 are given in the paragraph that follows Eq. ͑4͔͒. The calculation of the ͑angular͒ NACTs was carried out at the state-average completed active space selfconsistent-field ͑CASSCF͒ level employing the following basis functions: For the nitrogen we applied s, p, d, and f functions and for the hydrogens we employed s, p, and d functions all from the aug-cc-pVTZ set. We used the active space, including all seven valence electrons distributed on 10 orbitals. Six electronic states, including the three states specifically studied, were computed by the state-average CASSCF method with equal weights. In certain cases these calculations were repeated with seven states to check for convergence. The relevant ͑angular͒ NACTs, namely, 12 ͑q ͉ z͒ and 23 ͑q ͉ z͒ as well as the corresponding potential energy surfaces E j ͑q ͉ z͒; j =1,2,3 were calculated employing the MOLPRO program. 32 In what follows we discuss three adiabatic states of the HNH system, namely, the two ⌸ states that evolve to become 2 A 1 and 2 B 1 states and a third state which evolves from the ͑collinear͒ ⌺ state to become a second B state. In what follows we designate that second B state as B ⌺ . While doing that we show that as long as the region of interest is close enough to the axis of the molecule ͑for any z value and any geometry͒ the 2 A 1 and 2 B 1 states are strongly coupled to each other but only loosely coupled to other states of the molecule and therefore form, by themselves, a ͑two-state͒ Hilbert subspace. However, once this region is increased to include sections further away from the axis, these two states get coupled to other states and therefore fail to form a Hilbert subspace in the enlarged region. We intend to show that adding a third state, which in this case is the above mentioned B ⌺ state, leads to a Hilbert subspace formed by these three states. The main purpose of this article is to supply the numerical evidence for these assertions.

A. Potential energy curves
In Fig. 2 are presented the three potential energy curves under consideration as a function of q for the three different configurations. It turns out that the two ͑⌸͒ states 2 A 1 and 2 B 1 are the two lower ͑planar͒ states and that the 2 B 1 is the lowest of the two. 33 It is also well noticed that in all three configurations these two states are close to each other as long as q is rather small and eventually become degenerate once q = 0. As for the B ⌺ state, it is seen to be the highest among the three and that as long as q is small enough it is well separated from any state. As q increases it seems to approach the 2 A 1 state.

B. Nonadiabatic coupling terms
In Fig. 1 are presented the two major ͑angular͒ NACTs, namely, 12 which couples states 2 A 1 and 2 B 1 and 23 which couples the 2 A 1 state and the B ⌺ state. It is well noticed that along the HNH axis ͑where q ϳ 0͒ the values of 12 ϳ 1 but then as q increases the values of 12 decrease. The rate of decrease is uniform but z dependent. A different situation is encountered for 23 . Along the axis ͑where q ϳ 0͒ the values of 23 are very small ͑and, in fact, are expected to become zero, once q =0͒ but they increase as a function of q. The rate of increase depends on z and is not always uniform.
Of all three configurations the C 2v case ͑namely, when z = 1.95 a.u.͒ requires more attention, in particular, the second NACT, 23 ͑q ͉ z = 1.95 a.u.͒. As long as ͕0 Ͻ q Ͻ 1.1 a.u.͖ 23 is well defined as it couples the second state ͑that is, the state 2 A 1 ͒ and, energetically, the third state, which is the B ⌺ state. Unexpected things happen once q Ͼ 1.1 a.u. If again we assume 23 ͑q ͉ z = 1.95 a.u.͒ to be the NACT between 2 A 1 and the third ͑energetically͒ state we find that 23 ͑q Ͼ 1.1 a.u.͉ z = 1.95 a.u.͒ϵ0. The reason is that the third state, in that interval, is not the expected B ⌺ state but an A state. This, among other things explains why this NACT is identically zero. Continuing in the same manner we could establish that in along this interval the B ⌺ state is, energetically, only the fifth state. In other words it is 25 ͑q Ͼ 1.1 a.u.͉ z = 1.95 a.u.͒ which, along the interval q Ͼ 1.1 a.u., is the ͑natural͒ continuation of 23 ͑q ഛ 1.1 a.u.͉ z = 1.95 a.u.͒. More details on these and similar issues will be presented in a forthcoming publication.
The behavior of the RT NACTs is reminiscent on what is encountered in case of the JT NACTs. We found that at a given JT point of degeneracy all NACTs are zero except those that are formed at that particular point. 29͑a͒ Assuming the same to happen in the present case implies that since the ͑1,2͒ NACTs, 12 , are formed at points located along the HNH axis, the ͑2,3͒ NACTs, 23 , are expected to be zero along this axis. This simply means that the source for the 23 is located at other region͑s͒ in configuration space but then is forced to decay to zero while approaching the molecular axis.

C. The topological D-matrix and the single-valued diabatic potentials
In Figs. 3 and 4 are presented the diagonal elements of the D matrix calculated once for the two-state case ͓thus D ͑2͒ ͔ and once for the three-state case ͓thus D ͑3͒ ͔. The main issue in the present section is the comparison between the diagonal elements of these two matrices. In both cases we expect them ͑as a function of q but for various z values͒ to be 1. Whereas when they are close enough to the molecular collinear axis they are practically equal to 1, we start noticing deteriorations as q increases.
we show, in Fig. 3, only one of them. A different situation is encountered in case of D ͑3͒ . Here the three diagonal elements, D jj ͑3͒ ͑q ͉ z͒; j =1,2,3, are not necessarily identical and therefore the three of them are presented in Fig. 4. Comparing the curves in Figs. 3 and 4 one notices that the unit value for the various D-matrix elements is better preserved ͑as a function of q͒ in case of the three-state case.
A more direct comparison between the two-state case and the three-state case is given in Fig. 3 where a second curve, related to the three-state case, is presented. The values of this curve are defined as the geometric average of the three diagonal elements of D ͑3͒ , namely, It is obvious that for the cases that D jj ͑3͒ ͑q ͉ z͒ =1; j =1,2,3 we also have D ͑3͒ ͑q ͉ z͒ =1. From Fig. 3 it can be seen that as long as the calculations are done close enough to the molecular axis the D 11 ͑2͒ ͑q ͉ z͒ values are indeed close to 1. However, when the calculations are performed further away from the axis ͑i.e., when q becomes larger͒ D 11 ͑2͒ ͑q ͉ z͒ values deteriorate, sometimes even quite rapidly ͑except for the C 2v configuration where the deterioration is slower͒. A significant different situation is encountered for the D ͑3͒ ͑q ͉ z͒ values. They are not only close to 1 for small q values but remain 1 for a much larger q interval.
As mentioned earlier the diagonal elements of the D matrix contain information regarding the single valuedness of the diabatic potentials. Thus if the N ϫ N D matrix is diagonal then the corresponding N ϫ N diabatic potential matrix is single valued in the relevant region.
In general, we prefer N to be as small as possible. Therefore as long the region is close enough to the molecular axis the two-state ADT ͓see Eq. ͑5͔͒ is expected to yield a single valued 2 ϫ 2 diabatic potential matrix. 29͑d͒ However, as the region extends to larger q values, D ͑2͒ deviates significantly from the unit matrix while the three-state diabatization remains single valued ͑sometimes even for relative large q values͒ because D ͑3͒ is either equal to the unit matrix or very close to being a unit matrix. In other words if dynamic calculations have to be extended to such regions we have to employ the 3 ϫ 3 ADT matrix in order to form the required single-valued diabatic potentials. In the present study the three states are the two states 2 A 1 and 2 B 1 that evolve from the two degenerate ⌸ states and the B ⌺ state that evolves from the ⌺ state.

IV. ANALYSIS AND CONCLUSIONS
In this study we showed for the first time that the theory developed in recent years to study the topological effects formed by the JT degeneracy applies also in the case of the topological effects formed by the RT degeneracy. In order to achieve that, we had to define the molecular axis as the corresponding seam and refer to contours that surround this axis. In the present article we employed circular contours perpendicular to the molecular axis and centered at the axis; however, we expect the results to be unaffected if different contours are applied ͑as long as they surround the molecular axis, close enough to the axis͒.
Whereas the idea to use the ab initio treatment to calculate the terms that couple the two ⌸ states has been known for some time, 13 the way to incorporate the effect of the interaction due to other states seems to be new. The ability to incorporate the effect of RT coupling terms related to higher states opens up possibilities to predict the effect of indirect coupling terms due to remote states on processes that take place at lower states. An example for such a study is given here where we show how the interaction due to the "isolated" B ⌺ state ͑which evolves from the collinear ⌺ state͒ can be included.
One major issue has not been mentioned so far, namely, the possibility that the RT topological effects are affected by JT interactions. It is not obvious that these two interactions are connected because they are formed by two orthogonal angular motions, one is an in-plane ͑angular͒ motion and the other is an out-of-plane ͑angular͒ motion.
In Fig. 5 are shown the position of four degeneracy points responsible for the formation of the JT cis ͑recently we reported on such a study for the NHH configuration 34 ͒: two of them, namely, 12 JT ͑ , ͉ q , z͒ and 23 JT ͑ , ͉ q , z͒, that couple the two lower 2 A states ͑i.e., 2 A 1 and 2 A 2 ͒ and the two upper 2 A states ͑i.e., 2 A 2 and 2 A 3 ͒ are located along the C 2v line and the other two form two ci-twins of the type 23 JT ͑ , ͉ q , z͒ located on both sides of the C 2v line and couple the two upper 2 A states 35 ͑the two corresponding NACTs are of opposite signs͒. Here is the angular coordinate and is the radial coordinate of a point with respect to the position of the considered ͑JT͒ ci. Near each of the ci points is given the value of the corresponding topological ͑Berry͒ phase calculated according to the well known line integral which in case of a circular contour becomes 29͑c͒  It is noticed that in all four cases the ␣ phases are close to ± ͑and not 2͒ indicating that indeed the NACTs are formed by JT cis. A similar study, although for a spatial configuration determined by a different fixed R HH distances, supports the existence of four similar JT cis. 36 It is also noticed that the region that contains the four cis is the region where the RT interaction is relatively strong. Nevertheless the JT topological phases do not seem to be affected. FIG. 5. The position of the four JT cis in configuration space as found for R HH = 3.95 a.u.: ͑ࡗ͒ 12 JT ͑q ͉ z͒; ͑ϫ͒ 23 JT ͑q ͉ z͒. The C 2v ͑1,2͒ ci is located at a distance q = 2.27 a.u. from the collinear axis. The C 2v ͑2,3͒ ci is located at a distance q = 1.08 a.u. from the collinear axis. The two ͑2,3͒ ci twins are located at a distance q = 2.26 a.u. from the collinear axis and at a distance of 0.867 a.u. from the symmetry line. The numerical value near each ci stands for the topological phase ␣͑ ͉ q , z͒ as calculated for = 0.01 a.u. ͓see Eq. ͑20͔͒.
From the present study it is not clear yet whether the JT interaction interferes with the RT coupling ͑and vice versa͒ to the extent that it affects, for instance, the RT D-matrix elements. Although we have indications that these two types of interactions are transparent to one another still this issue has to be studied in more detail which will be done in forthcoming publications.
To summarize the findings of the present article we say the following: ͑1͒ The RT NACTs are similar to the JT ones; they are a result of a line of degeneracy points-seam-that yields a line of poles ͑called͒. ͑2͒ Their behavior and their effect on the diabatization process is well treated by the theory that originally was developed to study the JT NACTs.
Before concluding this article we would like to call attention to one important issue; while deriving the BO Schrödinger equation we anticipate it to contain both JT and RT singular coupling terms. Therefore the two have to be eliminated simultaneously in order to permit solving the resulting equations. This procedure, to the best of our knowledge, has never really been carried out.

ACKNOWLEDGMENTS
Three of the authors ͑M.B., Á.V., and G.J.H͒ acknowledge the US-Israel Bi-national Science Foundation for partly supporting this study. One of the authors ͑Á.V.͒ acknowledges the OTKA Grant Nos. T037994 and M041537 and the computational resources provided by the John-von-Neumann Institute, Research Centre Juelich ͑Project ID ehu01͒.