Generalized local frame transformation theory for ultralong-range Rydberg molecules

A detailed theoretical framework for highly excited Rydberg molecules is developed based on the generalized local frame transformation. Our approach avoids the use of pseudopotentials and yields analytical expressions for the body-frame reaction matrix. The latter is used to obtain the molecular potential energy curves, but equally it can be employed for photodissociation, photoionization, or other processes. To illustrate the reliability and accuracy of our treatment we consider the Rb$^*-$Rb Rydberg molecule and compare our treatment with state-of-the-art alternative approaches. As a second application, the present formalism is used to re-analyze the vibrational spectra of Sr$^*-$Sr molecules, providing additional physical insight into their properties and a comparison of our results with corresponding measurements.


I. INTRODUCTION
Rydberg molecules constitute one of the most exotic physical systems in quantum chemistry. They are effectively three-body systems composed of a ground state neutral atom (i.e. perturber), a positively charged core, and an electron; the subsystem formed by the latter two is a Rydberg atom. The Rydberg electron mediates an interaction that binds together the perturber and the ionic core. Early theoretical studies showed that the delicate nature of the binding mechanism results in the formation of a class of weakly bound molecules with bond lengths on the order of a few hundred nanometers [1,2]. These ultralong-range Rydberg molecules (ULRMs) are subdivided into "trilobite" [1] and "butterfly" [2] molecular species originating from the S− and P −wave "electronperturber" interactions. Although quite fragile, both UL-RMs were recently experimentally realized and observed [3][4][5]. One striking attribute of the ULRMs is that -despite the homonuclear nature of their constituents -they can possess huge dipole moments in the range of kilodebye [6,7]. Also, due to the resonant P −wave "electronatom" interaction, the butterfly molecules exhibit much deeper binding energies than the trilobite ones. The state of the art of these molecules is reviewed in Refs. [8][9][10].
From a theoretical viewpoint, the ULRM Hamiltonian possesses a fundamental attribute that has not received explicit attention: it exhibits local symmetries in different regions of configuration space. The exploitation of these local symmetries permits a more compact description of the physics of ULRMs and furthermore can be easily generalized to more complex scenarios, e.g. UL-RMs with multiple perturbers. One particular theoretical toolkit which can be employed in systems that possess * pgiannak@pks.mpg.de † meiles@pks.mpg.de ‡ robichf@purdue.edu § rost@pks.mpg.de local, rather than global, symmetries is the local frame transformation theory (LFT).
In the most general case these local symmetries can be associated with incompatible sets of approximately good quantum numbers which then are inter-related by the LFT. U. Fano in Ref. [11] and D. Harmin in Ref. [12] introduced the fundamental constituents of LFT theory in order to describe the Stark photoabsorption spectra of alkali Rydberg atoms. The LFT concept possesses a versatile scope capable of describing a plethora of different physical systems well beyond this original application. Its applications include, for example, studies of ultracold atoms [13][14][15][16][17], atoms in the presence of electric or magnetic fields [18][19][20][21][22] or generic trapping potentials [23], or Rydberg atoms with two valence electrons [24][25][26]. In the field of molecular physics, the LFT theory provides an insightful description of processes such as electronmolecule collisions [27][28][29], dissociative recombination of H + 3 [30,31], and Stark photoabsorption spectra in molecular Hydrogen [32]. In addition, the LFT theory has been extended to investigate the rovibrational spectra of diatomic molecules [33] or the electronic excitation of molecular ions [34]. However, for many years the LFT theory lacked a systematic pathway to improve its accuracy by including, for example, the physics of energetically strongly closed channels [35]. For this reason, the generalized local frame transformation theory (GLFT) was developed. It resolves the lack of closed channels in the original theory, providing a concrete framework based on more rigorous physical grounds.
In this study, we develop a non-perturbative theoretical framework for ULRMs which exploits the corresponding local symmetries. Our method combines the GLFT formalism with key ideas from multichannel quantum defect theory (MQDT) [36], going beyond previous studies [37,38]. MQDT permits us to treat the multichannel scattering of the Rydberg electron in the polarization potential of the neutral atom in a compact manner where all the relevant physics is encapsulated in the corresponding body-frame reaction matrix K. The GLFT treatment allows us to obtain K analytically without using the pseudopotential theory developed by Fermi and Omont [39,40]. In this manner, our theoretical framework avoids the limitations inherent to these methods [41,42]. For example, the GLFT approach guarantees the convergence of the potential energy curves within a particular electronic n Rydberg manifold without needing to add extra ones as in the case of the diagonalization method of Omont's pseudopotential theory.
We apply the GLFT method to the calculation of the potential energy curves of both Rb 2 and Sr 2 ULRMs. In the case of rubidium, we compare the GLFT theory results with those obtained by diagonalizing the Omont pseudopotential. Results from both methods for high electron Rydberg manifolds (n = 30) differ quantitatively, but agree qualitatively. However, at lower manifolds, i.e. n = 10, the diagonalization treatment exhibits even substantial qualitative differences from the GLFT approach. In the case of strontium, we investigate the ULRMs associated with the Rydberg s state since the vibrational spectra calculated using the GLFT can be compared directly to experiment [43,44]. This allows us to extract more precise information about the zero energy electron-Sr S−wave scattering lengthi.
The GLFT method gives identical results as the Green's function treatment developed in Refs. [2,45]. However, one major advantage of the GLFT is that its scope is much more general. The GLFT treatment encompasses all the relevant physics of the Rydberg atom interacting with the neutral one in terms of the bodyframe K−matrix which can be utilized to describe processes of predissociation, photoabsorption, or angular momentum changing collisions, as was pointed out also in Ref. [38]. Moreover, linking our method with the Coulomb-Stark frame-transformation in Ref. [35] permits the study of Rydberg molecules in external electric fields even in the regime where electronic n−manifolds are strongly mixed. Owing to the modular structure of the GLFT framework, the treatment presented here can also be easily extended to situations involving multiple perturbers providing a compact and accurate description of the corresponding potential energy surfaces.
This work is organized as follows: Section II lays out the concepts of our GLFT-MQDT theory for ULRMs. Section III is devoted to detailed comparisons with alternative approaches and ends with a reanalysis of Sr 2 experimental spectra. Since the derivation in Section II is quite involved, in our concluding section Section IV we provide a summary of the GLFT method and highlight the crucial expressions from the derivation before concluding.

II. THEORETICAL METHOD
A. The Hamiltonian of the three-body system and general considerations We are interested in the ULRM system, composed of a highly excited Rydberg atom A and a neutral ground state atom B and shown schematically in Fig. 1. To handle the multichannel scattering of the Rydberg electron by the ground state atom B and to obtain a body-frame K-matrix free from singularities, we will employ ideas from the GLFT and MQDT theory. The pair of quantum numbers l and m indicate the orbital angular momentum of the Rydberg electron and its projection relative to the Rydberg core A + , respectively. Similarly, the quantum numbers L and M denote the electron's orbital and azimuthal angular momentum relative to the perturber B, respectively. The internuclear axis ⃗ R = Rẑ is aligned with the z axis, and thus the azimuthal quantum numbers m and M are conserved, i.e. m = M . The position of the Rydberg electron is r with respect to the Rydberg core, and ξ = r − R with respect to the perturber (see Fig. 1). The full Hamiltonian in the frame of reference of the Rydberg core reads The term u sh (r) is a parameterization for the interaction between the Rydberg electron and the many electrons of the residual atomic core. This interaction is short-range, i.e. it vanishes for distances r > r A , and so its effect on the electronic motion is encapsulated by a set of quantum defects µ l . Note that in Eq. (1) and throughout this study we use atomic units unless otherwise specified. At distances r ≈ R, the Rydberg electron interacts with the perturber via a short-range potential,V B , of asymptotic form where α is the static polarizability of the neutral atom. At smaller ξ this interaction also includes electron correlation, interaction with the perturber's core electrons, and exchange interaction. As with the quantum defects, the effect of V B is imprinted on the wave function at large distance by a set of phase shifts. Since the interactions between the electron and the Rydberg core and perturber are limited to small volumes around each site, over the rest of the Rydberg's orbit the electron only experiences the Coulomb potential. In the spirit of MQDT, we postpone imposition of the physical boundary conditions to a later point in our derivation and impose "standing-wave" boundary conditions for now. This simplifies tremendously the analytic manipulations while keeping the scope of the treatment as general as possible. The Rydberg electron's wave function at r → ∞ reads where f c l (r) and g c l (r) are energy normalized regular and irregular Coulomb functions, respectively. Note that the symbolδ l ′ l denotes Kronecker's delta. The pair (f c l , g c l ) are evaluated at an electronic energy (R) which depends on the internuclear separation R, and are defined with the origin centered at the ionic core. K l ′ l (R) are the elements of the real, symmetric, body-frame K-matrix and satisfy Once the K-matrix is determined, we can impose the proper physical asymptotic boundary conditions by requiring that the electronic wave function vanishes at large r. This yields the determinantal equation where the principal quantum number ν = 1 −2 (R) is determined by the discrete electronic energy (R). These energies, parametrically depending on the internuclear separation, define the set of Born-Oppenheimer potential energy curves. The key quantity determining these potential curves is evidently the K-matrix. Therefore, the following subsections show how to analytically obtain the K-matrix using the GLFT.
B. Quantum-defect-shifted Coulomb functions and smooth K-matrix In the same spirit as Refs. [37,38] it is desirable to express the electronic wave function of Eq. (3) using an alternative pair of energy-normalized regular and irregular Coulomb functions, namely the so-called quantumdefect-shifted (QDS) Coulomb functions. This alternative pair of Coulomb solutions permits us to eliminate the potential term u sh (r) from Eq. (1) since its effects are encapsulated in the atomic quantum defects µ l . The energy-normalized QDS pair of solutions are related to the conventional regular and irregular Coulomb functions (f c l , g c l ) according to The QDS pair of solutions (F l , G l ) are associated with the K-matrix, not the K-matrix of Eq. (3). These matrices are inter-related via a simple matrix transformation (see discussion in Section II E). The electronic wavefunction from Eq. (3) at large distances is expressed as a linear combination of the QDS Coulomb functions, (7) and in general it fulfills the Lippmann-Schwinger equation whereĜ QDS C ≡ [ (R) + ∇ 2 r 2 + 1 r − u sh (r)] −1 is the quantum-defect-shifted principal-value Coulomb Green's function and V B is the electron-perturber interaction potential discussed above. This Green's function reads in terms of the QDS Coulomb functions as where r < = min(r, r ′ ) and r > = max(r, r ′ ). Note that G QDS C (r, r ′ ) obeys standing wave boundary conditions at large distances. The Lippmann-Schwinger equation in Eq. (8) provides us with the K-matrix, K l ′ l (R) = −π ⟨F l ′ (R) V B Ψ lm (R)⟩. As was shown in Refs. [23,35], by employing the Schwinger identity the K-matrix can be expressed in the following symmetric form: Similar to Refs. [23,35], the K matrix elements in Eq. (10) can be expanded over a complete basis set of wave functions φ LM ⟩ which are solutions of the electronperturber (e − B) subsystem. A key assumption, justified by the short-range nature of the e−B interactions, is that over the small region of space around the perturber where the polarization potentialV B prevails, the Coulomb field is approximately constant. The φ LM ⟩ wavefunctions fulfill the following Schrödinger equation: At distances ξ larger than the range of theV B potential, φ LM (ξ) acquires its asymptotic form where (f 0 LM (ξ), g 0 LM (ξ)) are the energy-normalized regular and irregular field-free solutions, namely the spherical Bessel and Neumann functions, respectively. The term δ L ≡ δ L (k) is the energy-dependent scattering phase shift induced by the polarization potentialV B . The relative momentum of the e − B subsystem is k ≡ k(R), and its R-dependence comes from the relation k(R) = 2 (R) + 2 R. The quantum numbers L and M label the angular and azimuthal momentum with respect to the perturber B and they are different from the quantum numbers l and m which are associated with wave functions centered at the ionic core A + . Also the azimuthal quantum numbers m and M are conserved; thus they remain the same for a wavefunction that is either centered at A + or B. The K-matrix in Eq. (10) expanded over the states Φ LM ⟩ takes the form Evidently, the computation of the K-matrix involves the evaluation of the matrix elements in Eqs. (13b) and (13c). However, the terms C T L ′ M,lm (R) and M L ′ L require the computation of complicated volume integrals which involve functions that are centered at two different locations of the configuration space: recall that the QDS function F lm (r; R) is centered around A + while the wave function φ LM (ξ) is centered around B. This challenge can be efficiently avoided by introducing the core idea of the GLFT approach, which enables us to transform the complicated volume integrals into simpler surface ones. In the following subsection we will introduce and derive the local frame transformation for Rydberg molecules, which inter-relates the functions F l ′ m (R)⟩ with the solutions f 0 LM ⟩. Such a transformation is valid only in a local region around B, where the motion of the Rydberg electron is mainly influenced by the polarization potential.

C. Local frame transformation for two different scattering centers
According to Fig. 1 the three-body system, A + − e − B, can be divided into two subsystems: (i) the Rydberg atom (A + − e), which is conveniently described by the QDS Coulomb functions centered on the Rydberg core, and (ii) the electron-perturber (e − B) complex that is addressed by free particle scattering solutions relative to the perturber B which is located at internuclear distances R from the Rydberg core. In the following, we derive the relevant expressions which connect the regular QDS Coulomb function F l ′ m (R)⟩ centered at A + with the energy normalized regular field-free solutions f 0 LM ⟩ centered at B.
As a first step, at distances r ≈ R we express the regular QDS Coulomb functions as a linear combination of the two linearly independent field-free regular and irregular functions that are centered at the Rydberg core A + : where the pair functions (f 0 lm , g 0 lm ) are the energynormalized spherical Bessel and Neumann functions respectively. The term W{⋅, ⋅} R indicates the Wronskian evaluated at R. It should be noted that Eq. (14) only applies in the vicinity around the perturber B.
The right-hand side of Eq. (14) is still expressed in terms of field-free regular and irregular functions relative to A + . These can be analytically re-expanded into fieldfree regular functions centered at a different point by exploiting basic attributes of bispherical harmonics [see Eq. (34) in Ref. [46]]. This identity relates the spherical functions, i.e. f 0 lm (r) and g 0 lm (r), to spherical Bessel where C j3m3 j1m1,j2m2 are the Clebsch-Gordan coefficients which account for the re-coupling of angular momenta relative to the Rydberg core with the corresponding angular quantum numbers that are associated with the perturber.
Substitution of Eqs. (15a) and (15b) into Eq. (14) results in the local frame transformation between the QDS regular Coulomb function and the field-free regular functions, Note that, when the quantum defects µ l vanish, Eq. (16) provides us with the transformation of the regular Coulomb functions f c lm (r) relative to A + in terms of the field-free regular solutions f 0 LM (ξ) relative to B. In addition, a similar local frame transformation can be derived for the irregular QDS Coulomb function defined in Eq. (6b) if needed. More specifically, as is discussed in the following subsection, the evaluation of the M L ′ L matrix elements is greatly simplified by frame transforming the Whittaker Coulomb functions which exponentially decay at large distances. The energy-normalized Whittaker Coulomb functions obey the relation w c lm (r) = Y lm (r)[g c l (r) + cot πνf c l (r)] and the corresponding local frame transformation reads: With this local frame transformation in hand, we compute in the following subsection the matrix elements necessary to determine the K-matrix.
D. Evaluation of the K-matrix in terms of the local frame transformation The corresponding volume integrals in the C T LM,lm (R) terms of Eq. (13b) contain the short-range polarization potentialV B . As was shown in Ref. [23], such integrals can generally be recast into surface ones. Given this, Eq. (13b) reads where dΩ ξ = sin θ ξ dθ ξ dϕ ξ is the solid angle differential element around the location of the perturber B. By employing the local frame transformation (see Eq. (16)), the C T LM,lm (R) matrix elements become The evaluation of the M L ′ L matrix elements in Eq. (13c) poses more challenges than the C T LM,lm (R) terms. The main source of difficulties arise from the QDS Coulomb Green's function which possesses a divergence for r → r ′ . However, in the spirit of the GLFT approach this type of divergence can be eliminated simply by adding and subtracting the field-free principalvalue Green's functionĜ 0 = [k 2 (R) 2 + ∇ 2 2] −1 which exhibits the same pole in the limit r → r ′ [35]. By forming the difference of principal-value Green's functions, i.e. ∆Ĝ ≡Ĝ 0 −Ĝ QDS C , the matrix elements M L ′ L are found to obey the following relation: In order to evaluate Eq. (20) it is desirable to isolate those terms in ∆Ĝ which depend on the atomic quantum defects µ l . This is achieved by expressing ∆Ĝ in terms of the Coulomb Green's functionĜ C which obeys exponentially decaying boundary conditions asymptotically, i.e. the physical Coulomb Green's function. This gives where All the functions that depend on the atomic quantum defects µ l are isolated in the term Λ, which is a finite sum which terminates at l N . For l > l N all quantum defects µ l vanish. Using Eq. (21) and the local frame transformation from Eq. (16), we recast the volume integrals into surface integrals as prescribed in Eq. (18) to obtain The terms From these closed-form expressions we can obtain the K-matrix free from any unphysical divergences. Since the physical boundary conditions are not imposed yet, the matrix elements of K are known for any angular momentum l. This is the focus of the following subsection, which addresses the importance of the energetically strongly closed channels, i.e. l > ν.
E. Pre-elimination of strongly closed channels in K-matrix and the connection relation to K−matrix In general, the motion of the Rydberg electron in the presence of a charged core and a neutral perturber can be viewed in terms of half collisions, as in the photoionization of Rydberg atoms. From this viewpoint, the (l, m) quantum numbers that are relative to the Rydberg core can be used to label the different collisional channels. The perturber couples only the different l-channels since m is a good quantum number, resulting in a non-diagonal K-matrix.
The different l−channels separate into two categories: (i) Weakly open channels, which refer to l-states with a weak centrifugal barrier yielding a classically allowed region between the Rydberg core and the perturber and (ii) Strongly closed channels, indicated by l−states associated with strong centrifugal forces imposing a classical forbidden region between the core and the perturber. In particular, we consider all l ≤ n * − 1 momenta in the weakly open channels, where n * is the integer lying in the interval ν < n * ≤ ν + 1. All angular momenta with l > n * − 1 are regarded as strongly closed channels.
For these strongly closed channels, the classically forbidden region causes the corresponding regular and irregular Coulomb functions to possess imaginary energynormalization constants; the electron's energy is below the combined centrifugal and Coulomb potential over the entire configuration space. Since we have not yet imposed the physical boundary conditions on these wave functions, the K-matrix contains those unphysical channels.
In order to properly treat the physics of the strongly closed channels, the wave function in Eq. (7) is partitioned into open ('o') and closed ('c') channels, as is usually done in multi-channel quantum defect theory. Sub-sequently, the closed channel part of the wave function is eliminated by imposing the physical boundary conditions asymptotically, i.e. by forcing the corresponding part of the wave function to exponentially decay at large distances. This pre-elimination of the strongly closed channels results in a wave function that involves only the weakly open channels and is associated with a real and symmetric physical K-matrix. In compact form the K-matrix reads: where the collective index o (c ) indicates all the angular momentum that fulfill the relation 0 ≤ l ≤ n * −1 (l > n * −1 ) for ν < n * ≤ ν + 1 where n * is an integer. In total there are N c closed channels. The matrix tan πν is diagonal with ν defined after Eq. (5). This elimination of strongly closed channels leads to the wave function As a final step, we compute the K-matrix using the K-matrix by linearly combining the solutions in Eq. (27) to construct the wave function that involves the regular and irregular Coulomb functions. This transformation can be expressed in a compact form: where cos πµ and sin πµ are diagonal matrices. Recall that µ l indicate the atomic quantum defects. Now that we have obtained a closed form expression for the physical K-matrix we can analyze the role of the strongly closed channels. Their importance is illustrated in Fig. 2. The Σ potential energy curves for Rb * −Rb near the n = 30 degenerate Rydberg manifold are shown for three different values of N c , N c = 10 (blue), N c = 20 (red), and N c = 60 (green). Note that the green dots represent the converged potential curve. Fig. 2 illustrates that at internuclear distances R < 1100 a 0 the trilobite potential curve is sufficiently converged by including a small number of strongly closed channels in the K phys matrix elements. However, for R ≥ 1100 a 0 , the convergence depends strongly on N c . Several numerical calculations for different electronic manifolds n show that the minimum number of strongly closed channels is N c ≈ 2n to ensure convergence of K phys out to and even beyond the classical turning point, R ≈ 2n 2 . In closing, we note that the theory of Refs. [37,38] neglects the strongly closed channels and therefore would not give accurate results for long-range Rydberg molecules.

A. Rubidium
In the following the GLFT approach is applied to the Rb * −Rb system to benchmark it against other state-of-the-art methods, namely the diagonalization of Omont/Fermi pseudopotentials and the Coulomb Green's function method. The S− and P −wave phase shifts of the "electron-perturber" subsystem are obtained by a non-relativistic two electron R-matrix method [47,48]. Recall that the S− and P −wave electronperturber phase shifts give rise to the trilobite and butterfly molecules, respectively. The atomic quantum defects and the "electron-perturber" phase shifts are the only auxiliary input parameters. Fig. 3 shows a comparison of the GLFT method (orange dots) and the diagonalization approach (blue solid line) where the electron-perturber interactions are modeled via Omont's pseudopotential [40]. For this paradigmatic calculation we consider the Σ potential energy curves relative to the n = 30 electronic Rydberg manifold. In the diagonalization approach the corresponding potential energy curves at n = 30 are calculated using one manifold above n = 30 and two additional manifolds below n = 30. On the qualitative level shown in Fig. 3(a), both methods are in reasonable agreement over the whole range. However, in the zoom-in plots in Fig. 3(b) and (c) near the trilobite and butterfly molec- ular curves, respectively, the quantitative differences between the GLFT and diagonalization approaches are evident. Note that the largest deviations are manifested near the avoided crossing between the trilobite and butterfly curve. Similar quantitative differences are seen in the quantum defect-shifted states (explored further in the following section) and, due to the extreme sensitivity of these weakly bound molecules, lead to noticeable disagreement between theoretical and experimental values for the binding energies. Note that the largest deviations are manifested near the avoided crossing between the trilobite and butterfly curve. Similar quantitative differences are seen in the quantum defect-shifted states (explored further in the Section III B) and, due to the extreme sensitivity of these weakly bound molecules, lead to noticeable disagreement between theoretical and experimental values for the binding energies. In the butterfly molecular curves shown in Fig. 3(c) the quantitative differences are more apparent, in particular, for internuclear distances to the left of the avoided crossing between the butterfly and the low-l potential curve "32p". The diagonalization method predicts that the wells in the butterfly curve are much shallower than those that are obtained within the GLFT framework. For example, at distances R = 500 a 0 the wells in the butterfly potential energy curve are ∼ 30% shallower than those of the GLFT approach whereas at shorter internuclear distances, i.e. R ≈ 100 − 200 a 0 , the corresponding deviation between the two methods increases at 46%. Such discrepancies in the relative depth of the wells strongly influences the number of vibrational states that they support. These differences between the diagonalization treatment and more sophisticated methods are known [41], and stem from the lack of a systematic pathway to increase the accuracy of the diagonalization method. The convergence of the potential energy curves at given n is not guaranteed by increasing the number of manifolds included in the basis.
This issue becomes more severe at low n. Fig. 4(a) and (b) illustrate the Σ potential energy curves of the n = 10 Rb * −Rb molecule. The differences between the GLFT approach (green dots) and the diagonalization method (blue and orange solid lines) are apparent even on a qualitative level. In this figure we illustrate also the convergence challenges inherent in the diagonalization method; the blue solid line in Fig. 4(a) and (b) is calculated by using the n = 8, 9, 10, 11 manifolds while the orange solid line is calculated using the n = 7 − 12 manifolds in the basis. These particular choices of the truncated basis yield potential energy curves that are closest to the corresponding GLFT calculations; nevertheless they are still quite different both from the GLFT calculation and from each other.
In addition, in panel (a) the diagonalization approach depends on the local momentum k(R) which is defined relative to the n = 10 manifold, i.e. k(R) = 2 R − 1 10 2 whereas in panel (b) the local momentum is defined by the relation k(R) = 2 R − 1 9.5 2 . This ambiguity in the local momentum is an important disadvantage, as the total electronic energy varies as a function of R and should be determined self-consistently (as is done in the Green's function method and in the present GLFT approach).
In Fig. 4 we observe that the butterfly potential curves obtained via the diagonalization method (see the blue and orange solid lines) are consistent with each other, however, they are inaccurate in predicting the range of internuclear distances R where the P −wave resonance induces a steep drop in the energy of the butterfly potential curve. For example, in Fig. 4(a) we observe that the butterfly potential plunges down at larger R for the diagonalization calculations than the corresponding ones in the GLFT method. By changing the definition of the local momentum k(R) in Fig. 4(b) the results of the diagonalization approach shifted to internuclear distances closer the GLFT calculations but the discrepancies between the two methods are still apparent. The robustness and the high level of accuracy of the present theory is illustrated in Fig. 5, where we compare the results of the GLFT theory (blue diamonds) with those from the standard Green's function method (orange dots) [2,45]. Note that Fig. 5 depicts both the Σ and the Π Rb 2 potential energy curves. The agreement between the two methods is excellent over the entire range of internuclear distances R up to the classical turning point, i.e. R = 2n 2 . Some deviations are evident beyond the classical turning point; these stem from the fact that the Green's function method which we used is not designed to correctly address the regime where the momenta k(R) becomes imaginary.

B. Strontium
Our discussion of the potential energy curves of Rb * +Rb Rydberg molecules emphasized their most prominent features, namely the potentials associated with the trilobite and butterfly states. These potentials, as a result of the high degeneracy of states in the Rydberg manifold, are typically many GHz deep. The remaining potential energy curves, associated with the non-degenerate quantum defect-shifted Rydberg states with low angular momentum, are the subject of the present section. The typical well depths are 10 − 100 MHz. The potential energy curve associated with a non-degenerate Rydberg state closely resembles the radial electronic wave function modulated by the S-wave and P -wave interactions with the perturber. Since we have confirmed in the previous section that the GLFT method circumvents the pitfalls inherent in attempts to increase the accuracy of the Fermi/Omont pseudopotentials and gives converged results, we now apply it to the calculation of the potential energy curves of Sr(ns)+Sr Rydberg molecules. These are ideal to study with the GLFT method because the commonly studied isotope, 84 Sr, does not have a hyperfine structure, the ns Rydberg states have no fine structure, and the spin-orbit splitting of the P -wave phase shift is expected to be small. As a result, the spin-independent theory developed here can be applied immediately, although we note in passing that the modular nature of the GLFT method can be readily extended to a spin-dependent Hamiltonian in future applications. By comparing the GLFT results with experimental signatures of these strontium Rydberg molecules, we can further benchmark the accuracy of the GLFT theory.
The vibrational spectra of these molecules are reported in Refs. [43,44] for principal quantum numbers ranging from n = 29 to n = 38. The n-dependence of the binding energies, E ∼ (n −μ s ) −6 , follows the anticipated Rydberg scaling law. Each potential curve supports around four vibrational states over this range of principle quantum numbers. These measurements form the basis for a series of experiments at increasingly high density and principal quantum number. In such regimes many perturbers can be located within the Rydberg orbit, and the experimental spectra exhibit features at energies which are sums of two or more dimer energies and are thus associated with polyatomic states of two, three, and more Sr perturbers [49,50]. At even higher densities and Rydberg levels the observed line shapes reveal details of the many-body response of the gas to the Rydberg impurity. Since the properties of a Rydberg atom embedded in these more exotic dense perturber scenarios are fundamentally linked to the energy levels of the strontium dimer, the GLFT calculations presented here serve a second purpose: to better calibrate these potential energy curves and the underlying atom-electron phase shifts.
As stated at the beginning of Section III A, to obtain potential curves for a new atomic species we only need to modify the quantum defects and scattering phase shifts which are the inputs to the GLFT calculation. The quantum defects of Sr are summarized in Ref. [8] and the electron-atom scattering phase shifts for Sr were calculated in Ref. [51]. It is clear from our GLFT calculations that these phase shifts are not accurate enough to obtain a theoretical spectrum in agreement with experiment. At large internuclear distances the calculated potential curves are too deep: the ground state molecules are more deeply bound in the calculation than observed in the experiment. At smaller internuclear distances the potential curves are insufficiently repulsive. The calculated excited states, particularly the 2-3 most highly excited ones, leak into the short-range region of strong ion-atom attraction, and their energies do not match the corresponding experimental measurements. It is not surprising that the calculated phase shifts at these very low scattering energies require tuning in order to obtain quantitative agreement (on the level of a few MHz at low n to even a few hundred kHz at higher n); this has been necessary in Rb as well. Rather than attempting a new fitting of the phase shifts, we adopt the same parameterization as in the theoretical analysis of Ref. [43]. There, the effective range theory for the energy-dependent S-wave scattering length a S (k) and a constant parameter for the energydependent P -wave scattering volume a 3 P (k) was used: The two fit parameters a S (0) and a P (0) are modified until the predicted vibrational spectra matches the measured spectra. Although this is a straightforward way to parametrize the energy-dependent phase shifts, we note that it has two key limitations. While a S (k) is the valid expansion of the scattering length in the k → 0 limit, as k increases a S (k) in Eq. (29) rapidly becomes only qualitatively accurate due to the presence of additional k 2 and k 2 ln k terms in the expansion. This means that although it can be used to obtain the zero-energy scattering length to a reasonable degree of confidence, it is not quantitatively reliable at smaller internuclear distances. This will be reflected in the binding energies of excited states probing these distances. On the other hand, a constant scattering volume is unphysical as k → 0, since the threshold law of the phase shift is quadratic in k as can be seen in Born approximation. This is typically not problematic in the calculation of Rydberg molecule potential curves since the P -wave contribution is much smaller than the S-wave interaction at distances R where this unphysical nature of a P is most prominent, i.e. in the vicinity of the classical turning point (R ≈ 2n 2 ). Although this parametrization could in principle be at least qualitatively accurate at larger k values, we note that the phase shifts of Ref. [51] exhibit significant energy dependence and this is only a crude approximation. We therefore emphasize that both of these parametrizations, but especially the effective Pwave scattering volume, are only convenient parametrizations for fitting and imply only very generic properties GLFT PT a S (0) a P (0) polarization potential (i) -12.65 9.6 (i) -13.   Table I. of the phase shifts except for the zero-energy scattering length, which can be fit quite accurately using only the ground vibrational state's binding energy. We obtain the phase shifts from Eq. (29) using Fig. 6 depicts two potential energy curves for the Sr(30s)+Sr Rydberg molecule. In the solid blue color we show the potential energy curve calculated using the GLFT approach with a S (0) = −12.65 a 0 and a P (0) = 9.6 a 0 . We refer to this model calculation as (i); within each model we calculate a family of potential curves for each n using these same parameters and computational technique. Importantly, (i) also includes the ion-atom polarization interaction, − α 2R 4 , where α ≈ 186 a.u. is the polarization of Sr. This attractive interaction is strongest at distances R < 500 a 0 and was neglected in Ref. [43]. The dashed brown curve in Fig. 6 shows the potential energy curve used in Ref. [43], which employed a model potential from first-order perturbation theory (PT) to Binding Energy (MHz) (b) 31 37 n  Table I treat S and P -wave electron-perturber scattering, (31) The semiclassical electron momentum is k s = 2 R − 1 (n −μ s ) 2 and Ψ s (R) is the Rydberg ns wave function. We refer to the model calculation using Eq. (31) with a S (0) = −13.2 a 0 and a P (0) = 8.4 a 0 as (ii). The a S (0) value in each case was obtained by fitting the binding energy of the ground vibrational state to the measured value. Since the P -wave interaction and ion-atom interaction are very weak at such large R values this state's energy is essentially fixed by the strength of the S-wave interaction. After fitting a S (0), a P (0) was fit by aligning the theoretical and experimental binding energies for the higher excited states. The vibrational wave functions bound in potential (i) are shown also in Fig. 6 in yellow. Those obtained using potential (ii), not shown, are nearly identical and have very similar binding energies, despite the differences in the potential energy curves. Since the P -wave interaction in both models is repulsive enough to restrict nearly all of the vibrational states to distances larger than 500 a 0 , the polarization potential has only a small effect. Fig. 7(a) compares the binding energies computed using the GLFT approach (model (i)) with the measurements reported in Ref. [43]. This demonstrates that the fitted a S (0) and a P (0) values used in model (i) are sufficient to closely reproduce the measured spectra over this whole range of n values. To verify that the fit parameters are independent of n we chose to fit them to the vibrational energies using only the 30s spectrum. The filled orange circles mark states whose line strengths are roughly equal, while the open purple circles mark states whose line strengths are less than 20% as strong as the average orange point. The line strengths are computed using a simple approximation for the Franck-Condon factor assuming a flat initial scattering state between the two atoms, as in Ref. [43]. Fig. 7(b) shows the results of vibrational states using the different models which are summarized in Table I. More specifically, in Fig. 7(b) we compare only the 30s data, but now also for a perturbative calculation using Eq. (31) including the polarization potential. We label this calculation as model (iv). By this additional comparison we are able to better understand the way in which the fit parameters compensate for inaccuracies in the 30s potential curves, whether these stem from the choice of method (GLFT or PT) or from the inclusion of the polarization potential. To obtain a theoretical spectrum with model (iv) matching the experimental data when including the polarization potential, we had to re-fit the a P (0) parameter to a P (0) = 9.8 a 0 . Although the model calculation (ii) used in Ref. [43] apparently ignored the polarization potential, this comparison shows that this error was compensated for by the difference ∆a P = 1.4 a 0 between models (ii) and model (iv). Fig. 7(b) shows also the results of model (v) which refers to a GLFT calculation using the scattering parameters of model (ii) and neglecting the polarization potential. This method clearly fails to predict the measured binding energies.

Lessons learned from different models
From these calculations, we arrive at three conclusions. First, that the GLFT model using fitted effective electron-atom scattering phase shifts (model (i)) reproduces, for all n values where measurements exist, the experimental spectra. Complementing the theoretical comparisons in Section III A, this provides an experimental confirmation of the validity of the GLFT approach. Some minor discrepancies are, however, visible in the comparison. The theory predicts weakly bound vibrational states (binding energies of 3MHz or less) with weak but finite line strength which are not observed in the experiment. These might be suppressed in the experiment by additional contributions to the line strength not accounted for in our approximate Franck-Condon factors, or they could be obscured by the atomic resonance, or indeed they could signal an additional energy dependence that cannot be compensated for by fitting a constant P -wave scattering volume or which stems from the higher-order effective range terms neglected in a S (k). A more glaring discrepancy is that the theory predicts two deep bound states for the 29s state with strong line strengths which are not observed. The theory calculations of Ref. [43] using model (ii) show also these extra vibrational bound states; based on the trends of the other vibrational states as a function of n one would expect these states to exist. Further study of the experimental spectrum is necessary to resolve this issue and identify if it has an origin in the experimental setup or if it heralds additional physics not included in the theory.
Our second conclusion is that the electron-strontium scattering length must be in the neighborhood of −12.65 a 0 . We estimate an uncertainty of 0.1 atomic units on this value based on the fit of a S (0). The extracted scattering length from model (ii), -13.2 a 0 , is about 5% different; this difference can be attributed to the known differences between calculations using the Fermi pseudopotential and the GLFT calculation. Indeed, in the trilobite potentials dominated by S-wave scattering presented in our study of Rb (see Fig. 3(b)) differences on this order were already visible. The a S (0) values obtained using the models (i) and (ii) both differ quite strongly, by 25%, from the calculated scattering length [51]. Although the authors of Ref. [43] attributed this to the uncertainty in the potential energy calculation, we can now claim on more rigorous theoretical grounds that the greater part of this difference is due to a real discrepancy between the scattering length provided in Ref. [51], and only the 5% error between −12.65 a 0 and -13.2 a 0 is due to the approximations made in the perturbative potential energy curve of model (i). That the calculated scattering length is overestimated is not surprising due to the challenges of converging the atom-electron scattering calculations at such low energies.
Third, we conclude that the P -wave scattering volume must be large and positive in order to produce a repulsive barrier at short-range which suppresses the effect of the polarization potential and localizes the vibrational states in the outer wells. This implies that the P -wave phase shift must be negative and significantly larger in magnitude than given by the calculations, particularly at high k values. These results show that further work is necessary in order to produce a more accurate set of scattering phase shifts. Caution is also warranted since, as seen in Figs. 6 and 7, the scattering parameters and even the potential energy curves themselves are not uniquely determined by the binding energies, and fitting of the binding energies alone is likely insufficient without additional theoretical calculations of the strontium-electron phase shifts to pin down these phase shifts.

IV. SUMMARY AND CONCLUSIONS
A. Summary of the method Since the derivation of the GLFT approach for Rydberg molecules involves many steps, we provide here a summary of the key equations that must be implemented in order to utilize this method. Also, note that Fig. 8 outlines the summary of these key expressions. First, Eqs. (13a) to (13c) define the K matrix which is associated with the QDS energy-normalized regular and irregular Coulomb functions. Using the Schwinger identity, the K matrix acquires a separable form and depends on two terms: the C T LM,lm (R) and the M L ′ L matrix elements. These are expressed in terms of the local frame transformation in Eq. (19) and in Eqs. (23) to (25), respectively. Eq. (16) provides the local frame transformation relation which connects the QDS Coulomb functions centered at the Rydberg core with the regular field-free functions around the perturber.
In the spirit of MQDT theory, the strongly closed channels of the K matrix are eliminated. This elimination step yields the physical K matrix, namely K phys , which is given in Eq. (26). In a final step, the physical K phys is connected to the K-matrix via the expression Eq. (28). Note that K−matrix is associated with the energy-normalized pair of regular and irregular Coulomb functions. For a specific Rydberg molecule two inputs are required: the atomic quantum defects µ l and the electronatom scattering phase shifts δ L . Using these two parameters, the K-matrix is inserted in the determinantal equation Eq. (5). The roots ν(R) of Eq. (5) determine the molecular energy curves via (R) = − 1 2ν 2 .

B. Conclusions
Based on the generalized local frame transformation theory we have developed a formalism to describe asym-metrically excited three-body systems, exemplified here with an excited Rydberg atom interacting with a ground state atom. By employing the key concept of the local frame transformation, we have obtained closed form analytical formulas for the body-frame K−matrix associated with a diatomic ultra-long-range Rydberg molecule. The potential energy curves are then obtained by numerically solving the one-dimensional determinantal equation expressed in terms of the K-matrix. We have shown that this GLFT approach provides potential energy curves which are known to be more accurate than those obtained via the diagonalization approach, which cannot be rigorously converged. We have used this advantage of the GLFT to re-analyze the vibrational spectra of strontium Rydberg molecules, and have seen that it can be used to extract a more accurate zero energy scattering length from experimental measurements. One major advantage of the GLFT method over the Green's function treatment, which operates at the same level of accuracy, is that it can be easily extended to other physical systems, for example to Rydberg atoms in the presence of multiple perturbers. This would provide quantitative improvements to the theory of polyatomic Rydberg molecules and Rydberg composites developed in [50,52,53]. In addition, due the modularity of the GLFT toolkit, the local frame transformation for Rydberg molecules presented here can be combined other frame transformations in order to investigate more complicated physical systems, e.g. Rydberg molecules in external electric fields or include relativistic effects.