Multichannel quantum-defect theory for ultracold atom-ion collisions

We develop an analytical model for ultracold atom-ion collisions using the multichannel quantum-defect formalism. The model is based on the analytical solutions of the r^-4 long-range potential and on the application of a frame transformation between asymptotic and molecular bases. This approach allows the description of the atom-ion interaction in the ultracold domain in terms of three parameters only: the singlet and triplet scattering lengths, assumed to be independent of the relative motion angular momentum, and the lead dispersion coefficient of the asymptotic potential. We also introduce corrections to the scattering lengths that improve the accuracy of our quantum-defect model for higher order partial waves, a particularly important result for an accurate description of shape and Feshbach resonances at finite temperature. The theory is applied to the system composed of a 40Ca+ ion and a Na atom, and compared to numerical coupled-channel calculations carried out using ab initio potentials. For this particular system, we investigate the spectrum of bound states, the rate of charge-transfer processes, and the collision rates in the presence of magnetic Feshbach resonances at zero and finite temperature.


Introduction
The successful realization of systems combining ultracold atoms and ions is stimulating an increasing interest in the physics of atom-ion collisions at very low collision energy. The first experiments were realized for clouds of atoms and ions stored in hybrid dual charged-neutral traps at mK temperatures [1,2]. More recent experiments performed with atoms in the micro and nano Kelvin temperature domain consist in immersing single ions in a Bose-Einstein condensate [3][4][5]. Understanding atom-ion collisions in the quantum regime is an essential elementary step before a more complete many-body description of these systems can be developed. While atom-ion collision properties are well known for high collision energies [6,7], a theoretical description in the ultracold domain is still under development [8][9][10][11]. Atom-ion and neutral atoms scattering are significantly different, mainly due to the relatively long-range character of the polarization interaction between atom and ion. Beyond the fundamental interest of their collision properties, these systems are also very attractive for quantum information processing [12,13]. Hybrid architectures may profit from advantages offered by both ionic and atomic species, namely a short computation time for charged particles and a long-coherence time for neutral atoms. Ultracold charged-neutral systems are also expected to exhibit interesting many-body effects, including for instance nontrivial modifications of the condensate wave function in the presence of ionic impurities [14] or the creation of mesoscopic size molecular ions [15].
In a previous paper [10] we have studied the basic properties of atom-ion scattering and bound states, using 40 Ca + − 23 Na [16] as a reference system. We have developed an effective atom-ion collision model by applying the multichannel quantum defect theory (MQDT) [17][18][19] to the polarization potential, which scales as r −4 at large atom-ion distance r. We have verified that MQDT model predicts very accurately all atom-ion scattering properties at ultracold temperatures by comparing to numerical close-coupled calculations performed on ab initio 40 Ca + − 23 Na potentials. In the literature, the MQDT approach has already been applied to describe scattering and bound states in electron-ion core [17], electron-atom [20] and neutral atom systems [21]. Our model of atom-ion collisions combines knowledge of the analytical solutions for the r −4 asymptotic potential [20,22,23] with the idea of a frame transformation [21,[24][25][26][27][28]. The latter, applied at small distances allows a reduction of the number of quantum-defect parameters needed to represent the effect of the short-range interaction potential.
In [10] we have applied MQDT to the study of magnetic Feshbach resonances and of radiative charge exchange processes, including the effects of Feshbach and shape resonances. In this paper we present a detailed derivation of our atom-ion model, based on the Mies formulation of MQDT [19]. We begin by introducing the quantum-defect approach, and discuss the properties of analytical solutions of the r −4 asymptotic potential. We derive asymptotic expansions of the model quantum-defect functions, valid for energies below a characteristic quantity E * defined in terms of the atom-ion interaction strength. For energies above E * , an efficient numerical algorithm to determine the characteristic exponent and the solutions of the Mathieu equation needed for the MQDT is proposed. Knowledge of the quantum defect parameters is crucial for understanding atom-ion collisions in the ultracold domain, and allows us to derive several non trivial analytical results concerning the behavior of weakly bound states, s-wave magnetic Feshbach resonances, and radiative charge-transfer probabilities. In comparison to the model presented in Ref. [10] we introduce an angular-momentum-dependent correction that improves the accuracy of MQDT for high order partial waves. This modification allows the regime of applicability of our model to be extended up to mK temperatures. Using this revised model and numerical close-coupled calculations, the effect of finite energy Feshbach resonances on elastic and charge-transfer collision rates is investigated as a function of the magnetic field.
In the present paper we also include details on the semiclassical description of radiative chargetransfer rates. A combination of the semiclassical approximation and MQDT scaling functions allows us to represent the charge-exchange loss rate as the product between the classical Langevin result and a universal quantum correction dependent only on a single parameter, the short-range phase. In [10] we discovered that the main product of a charge-transfer process for 40 Ca + − 23 Na are molecular ions. Here we present a more detailed analysis by analyzing the population of vibrational states formed after the electron transfer. A semiclassical formula describing this distribution is derived. The very good accuracy of our MQDT analytical results is tested by comparison with the numerical solution of the close-coupled Schrödinger equation for the 40 Ca + − 23 Na system.
The paper is organized as follows. Section II introduces the MQDT formalism. Section III discusses the properties of scattering and bound states for a single channel potential, based on the analytical solutions for the polarization potential. In particular, section III.A derives the analytical solutions, section III.B presents the results of the scattering states and analyzes the applicability of the semiclassical description. Section III.C presents the analytical results for weakly bound states and section Section III.D studies the low-energy behavior of the MQDT functions. The idea of the frame transformation is discussed in Section IV, using a particular example of an alkali atom and an alkali-earth ion. Section IV.A generalizes the frame transformation to the presence of magnetic field, and Section IV.B derives the corrections to the quantum-defect matrix for higher partial waves. Section V describes the radiative charge transfer: the perturbative approach based on the Fermi golden rule, and the semiclassical description of the charge-transfer probability at short range. The theoretical model is applied to the 40 Ca + − 23 Na system in Section VI. Section VII presents the conclusions, and five appendices give some technical details related to the derivation of the analytical solutions and the multi-channel calculations of the radiative charge transfer.

Multichannel formalism and quantum defect theory
In this section we summarize the basic formalism of the multichannel quantum defect theory, adopting the formulation introduced by Mies in Ref. [19]. Atom-ion collisions are described by a N-channel close-coupled radial Schrödinger equation Here µ = m i m a /(m i + m a ) denotes the reduced mass, and W(r) and F(r) are N × N matrices representing the interaction and the radial solutions, respectively. A general solution to the Nchannel scattering problem is given by a set of N linearly independent wave functions where |Φ j are channels states describing the internal spin degrees of freedom, Y ℓ j (r) denotes the angular part of the solution (spherical harmonic) for the channel j. The interaction matrix is asymptotically diagonal where E ∞ i are the threshold energies for the molecular dissociation, ℓ i is the partial wave quantum number of channel i, C 4 = αe 2 /2 with α denoting the static dipolar polarizability of the atom and e is the ion charge. Here we neglect the contribution of the higher-order dispersion terms to the long-range potential, which give relative small corrections in atom-ion scattering [10].
Given the total energy E, the channel states can be classified as open for E > E ∞ i or closed for E < E ∞ i . In the former case the asymptotic wave number k i = 2µ(E − E ∞ i )/ 2 is real and positive, while in the latter is purely imaginary k i = e iπ/2 |k i |. The solution matrix F(r) can be split into blocks Here, J ij (r) → δ ij sin(k i r − ℓ i π/2)/ √ k i and N ij (r) → −δ ij cos(k i r − ℓ i π/2)/ √ k i (r → ∞) exhibit asymptotic behavior associated with the spherical Bessel functions j ℓ (kr) and n ℓ (kr). The constant matrix A oo depends on the boundary conditions at r → ∞. In particular, the choice A oo = (1 − iK oo ) −1 corresponds to the usual incoming-wave boundary conditions, where the amplitude of the outgoing wave is determined by the scattering matrix S oo = (1+iK oo )(1−iK oo ) −1 .
The basic idea of the quantum-defect theory is to introduce a set of parameters describing the short-range behavior of the wave function, that weakly depend on total energy E, and can be used to predict the system properties as E crosses the dissociation thresholds of the individual channels. As we will show later, in the ultracold domain the quantum defect parameters for atomion collisions are also weakly dependent on the relative orbital angular momentum ℓ and can be taken as constant, at least for the lowest partial waves.
We now specialize the formulation of MQDT developed by Mies for atomic collisions [19] to our atom-ion system. The starting point is the choice of a set of reference potentials {V j (r)} that should reproduce the asymptotic behavior of the interaction matrix at large distances V j (r) r→∞ −→ W jj (r) but can otherwise be arbitrary. One associates to the reference potentials V i (r) a pair of linearly independent solutionsf i (r) andĝ i (r) that have short-range WKB-like normalization The amplitude α i (r) fulfills the inhomogeneous Milne equation: [d 2 /dr 2 +k i (r) 2 ]α i (r) = α −3 i (r) [29], with the local wavevector k i (r) = 2µ(E − V i (R)/ , while the phase dβ i /dr = 1/α 2 i . Since the reference potentials {V j (r)} reproduce the asymptotic behavior of W(r), the exact solution to Eq.(1) can be expressed at large distances in terms of pair of functionsf(r) ≡ {δ ijfi (r)} and Here, Y is the quantum-defect matrix that play a central role in the MQDT analysis. In contrast to the S and K scattering matrices, Y remains analytic across the thresholds, and has in general only a weak dependence on energy.
The observable properties depend on the asymptotic behavior of the solutions at large distance. It is therefore convenient to introduce another pair of solutions, f i (r) and g i (r), together with the physically well-behaved solution φ i (r) for the closed channels, that have energy-like normalization as r → ∞ The factor 1 2 in the last equation is introduced for convenience in order to simplify the value of Wronskian [19]. Here, ξ i is the scattering phase shift for a channel i. The MQDT functions C i (E), tan λ i (E) and ν i (E) relate the solutions (8) to (6), and are defined as follows Now, when all the channels are open: Using the relationships (9a)-(9c) one can show that where C(E) ≡ {δ ij C i (E)} and λ(E) ≡ {δ ij λ i (E)}. At high energies the WKB approximation is valid at all distances, and the functionsf ,ĝ become identical to f and g, respectively. Therefore, in this situation the MQDT parameters behave like C i (E) → 1 and tan λ i (E) → 0, and Eq. (11) reduces to R(E) ∼ = Y(E). Next, with the help of Eq. (5) one can relate the reactance matrix K and the scattering matrix S to the matrix R The current derivation has to be modified in the presence of closed channels: E < E ∞ i for i = N o + 1, . . . , N. In this case, we impose the requirement that the wave function of the closed channels decays exponentially for large r, i.e. the closed channel wave functions are proportional to φ i (r). This results in the renormalization of the open-open block of the quantum-defect matrix [19] where ν(E) ≡ {δ ij ν i (E)} The scattering matrices can now be calculated from Eqs.

Analytical solutions
Here we focus on the Schrödinger equation for a single channel where we include only the long-range part of the atom-ion interaction and the centrifugal barrier for partial wave ℓ In the following we will work in dimensionless units, where the length is expressed in units of R * ≡ 2C 4 µ/ 2 and energy in units of Table 1 presents the characteristic lengths R * and energies E * for some sample combination of alkali atoms and alkali-earth ions. Fig. 1 shows the long-range atom-ion potentials for the lowest partial waves, where the squares mark the top of centrifugal barriers occurring at r max = √ 2/ ℓ(ℓ + 1)R * with energy E max = 1 4 ℓ 2 (ℓ + 1) 2 E * . The Eq. (16) in characteristic units of R * and E * takes the form This equation can be solved analytically [22,23,30] by substituting F (r) = ψ(r)r 1/2 and r = e z E −1/4 , which yields the Mathieu's equation of imaginary argument with a = (ℓ + 1 2 ) 2 and q = √ E. Although the Mathieu's equation is well known in mathematical physics (see e.g. [31,32]), we include a brief derivation of its basic properties in the Appendix Appendix A. The Appendix also discusses the basic methods we applied for its numerical solution.
We denote two linearly independent solutions to Eq. (17) by T ν (r) and T −ν (r), where ν is the complex number called the characteristic exponent (see Appendix Appendix A for definition). In view of the asymptotic properties of T ν (z) discussed in Appendix Appendix A), we can construct two solutionsf (r) andĝ(r) defined by Eqs. (6) where and φ is some parameter that can be interpreted as the short-range phase. As can be easily verified, f (r) andĝ(r) exhibit at distances r ≪ 1 semiclassical behavior given by Eqs. (6) with α(r) ∼ = r and β(r) ∼ = −1/rf The asymptotic behavior given by Eqs. (22) can be also obtained by solving directly Eq. (17) with centrifugal and energy terms neglected. Since r sin(−1/r + φ) is an exact solution of (16) for ℓ = 0 and E = 0, it can be used to derive a simple relation between the short-range phase φ and the s-wave scattering length a = R * cot φ, that follows from the asymptotic behavior: F (r) → r − a (r → ∞) at E = 0. At large distances, the solutionsf (r) andĝ(r) behaves according tô where Here, η = π 2 (ν − ℓ − 1 2 ), m ν = S ν (4/q) ν , and S ν is a function of ν, that is defined in terms of a continued fraction (see Appendix Appendix A for the definition).
Making use of the asymptotic behavior (24) one can find the phase shift ξ and all MQDT functions defined in (9a)-(9c) where χ = √ −E, andξ is the phase shift of a second short-range normalized solutionĝ(r) Finally we point out that formulas (24) describing the asymptotic behavior off (r) andĝ(r) can be used to determine the scattering matrix K directly from the quantum-defect matrix Y where We note that in principle the parametrization in terms of the quantum-defect matrix Y(φ) depends on the short-range phase φ that is fixed by the choice of the reference potentials. In contrast, the scattering matrix K will depend only on the actual scattering lengths associated with the scattering channels and should be independent of a particular set of reference functions that determine specific φ for that choice.

Scattering states and semiclassical approximation
Using Eq.(26) and the small-q expansions presented in Appendix Appendix A.4 one can derive the well known threshold behavior of the phase shifts in the polarization potential [23,30] tan Here, we have used the relation (23) between the short-range phase ϕ and the scattering length a of the reference potentials, in order to express all the formulas in terms of a. Because of the long-range character of the polarization potential, the scattering length can be defined only for ℓ = 0, while for ℓ > 0 the threshold behavior of the phase shifts is dominated by the contribution of the polarization potential. Fig. 2 shows the energy-dependent scattering length a(k) = − tan ξ 0 /k as a function of φ for diferent values of the energy. The case k = 0 corresponds to a = R * cot φ.
At sufficiently high energies the scattering can be described using the semiclassical approximation. To identify the crossover from the quantum to the semiclassical regime, we plot in Fig. 3 the quantity λ ′ (r)/2π, where λ(r) is the local de Broglie wavelength. A necessary condition for the applicability of the WKB approximation, expressed in terms of the local wavelength, is λ ′ (r) ≪ 2π. We observe that this condition is first violated at distances comparable to R * , while the WKB approximation remains valid at small and large distances. We have verified that to a good approximation the wave function can be calculated within the semiclassical approximation when λ ′ (r)/2π 1/2, a condition fulfilled for energies E 25E * . Figs. 4 and 5 compare predictions of the WKB method with the exact wave functionsf (r) for two values of energy: E = E * and E = 100E * . In addition they show the small-r and large-r asymptotic formulas (22a) and (24a), respectively, and derivative of the local wavelength: λ ′ (r)/2π.

Bound states
For negative energies the exponentially decaying solution is given by a linear combination off (r) andĝ(r) [see Eq. (9c)] with ν given by (29). In the quantum-defect approach the bound state spectrum is determined by the condition Y + tan ν(E) = 0, which is a single-channel version of Eq. (15). In the presence of a single channel one can take Y = 0, which is satisfied when the reference potential is chosen to reproduce the same scattering length as the real physical potential. This yields the condition tan ν(E) = 0. In this case the bound state wave function φ(r) of Eq. (9c) contains only thef (r) component, with the short-range phase φ given in terms of the scattering length by Eq. (23). The condition determining the energy of a bound state can be rewritten as with S ν (E) defined in the Appendix Appendix A.3. Fig.6 shows energies of the bound states versus the short-range phase φ for the lowest partial waves. We observe that for ℓ even the bound states disappear at the threshold at φ = 0 (a = ±∞), while for odd ℓ this happens for φ = π/2 (a = 0). Thus, for a polarization potential there is an ℓ = 2 periodicity in the values of φ at the threshold, which is reminiscent of ℓ = 4 periodicity for van der Waals potential [33]. The energies of the bound states for a = ±∞ and for different partial waves are shown in Fig. (7). For even ℓ this defines characteristic energy bins, which determine the positions of the last bound states in the spectrum. For instance, the last s-wave bound state is located in the energy range E/E * = [−106, 0] Applying small q expansions of ν and S ν (see Appendix A) one can expand the right-hand-side of Eq. (34), which yields for s-waves and for p-waves  with κ = √ −E. Inverting the former equation, we can write the energy of s-wave bound states in powers of 1/a The first term on the right-hand-side is the universal energy of a weakly bound state, the higher order term represents the correction that is specific for 1/r 4 potential. Fig. 8 shows the binding energy for an s-wave bound state versus the inverse of the scattering length. We observe that the range of the applicability of the universal formula E = −1/a 2 is very narrow. Inclusion of the higher order correction in 1/a given by Eq. (37) improves the agreement with the exact one. In fact the approximation (37) works reasonably for |E| 0.5E * .  of φ in Fig. 2. The calculations are performed under the assumption that the short range phase φ does not depend on ℓ, such that all the MQDT functions can be parameterized with a single ℓ-independent φ. Substituting the small-q expansions given in Appendix (Appendix A.4) into Eqs. (27) and (28) one obtains

Behavior of MQDT functions for polarization potential
Here, δ ℓ,0 = 1 for ℓ = 0 and is zero otherwise. At sufficiently large energy, the semiclassical solution begins to be valid at all distances, implying C(E) → 1 and tan λ(E) → 0. One can verify that the MQDT functions presented in Figs. 9-11 exhibit the small-q and large-q asymptotic behavior given by Eqs. (38). For nonzero angular momenta one can observe the existence of sharp peaks in the function C(E), corresponding to shape resonances, that appear due to the presence of quasi-bound states behind the centrifugal barrier.

Frame transformation
In the ultracold regime variations of the total energy E are typically much smaller than the depth of the potential at short range, where the matrix Y is defined. Therefore, it is justified to neglect  the dependence of Y(E) on energy and to set Y(E) ∼ = Y. In this way, the matrix Y can be determined at a single value of energy, for instance from the reaction or scattering matrices K or S by applying formulas (12) and (13), respectively. An alternative way is to apply the frame transformation technique [26,34] that provides a very convenient way of parameterizing the short range matrix Y in terms of few parameters, e.g. scattering lengths.
In this paper we focus on the collisions of an alkali atom with an alkali-earth ion in their electronic ground states. Hence, the asymptotic channel states can be characterized by the hyperfine quantum numbers f 1 ,m f 1 and f 2 ,m f 2 for ion and atom, respectively, and by the angularmomentum quantum numbers ℓ and m ℓ of the relative motion of the atom and ion center of masses. We label the internal quantum numbers in the asymptotic basis by α = {f 1 f 2 m f 1 m f 2 }, hence the asymptotic channels are characterized by i = |αℓm ℓ . At short distances it is more convenient to characterize the channel states in terms of the total electron spin S = s 1 + s 2 and total nuclear spin I = i 1 + i 2 , where s 1 , s 2 are electron spin of ion and atom respectively, and i 1 , i 2 denote their nuclear spins, respectively. This basis is characterized by two additional quantum numbers, the total hyperfine angular momentum F = f 1 + f 2 = I + S and its projection M f on the axis of quantization. We denote the internal quantum numbers in the molecular basis by β = {ISF M F }, so the molecular channel states are characterized by j = |βℓm ℓ . The frame transformation is a unitary transformation between channels α and β, that can be written as Here, for simplicity we have omitted the quantum numbers ℓ and m ℓ that are conserved if one neglects weak dipolar interactions. Symbol (m f 1 m f 2 |F M F ) stands for the usual Clebsch-Gordan coefficient while is the transformation between f 1 f 2 and IS coupling schemes, given in terms of the Wigner 9j symbol.
The knowledge of the analytical solutions for the polarization potential suggest that the reference potentials can be chosen to contain only the diagonal long-range part of the interaction matrix: According to Eqs. (22a) and (22b) the wave functions in pure polarization potential are singular at r → 0. In this case the standard boundary conditions F(r) → 0 (r → 0) that is imposed on the physically meaningful solution, has to be replaced by the boundary condition (22a) with some short range phase φ i .

Angular-momentum-insensitive quantum-defect matrix
The exchange interaction that mixes asymptotic channel states α takes place usually in some range of distances R 0 of the order of few tens of a 0 . At such distances the interaction potential between species is much larger than the hyperfine or Zeeman splittings in the presence of an external magnetic field. Therefore it is convenient to define the matrix Y at distances r R 0 , where the exchange interaction is negligible, and W ij (r) ∼ = δ ij C 4 /r 4 for R 0 r ≪ R * . This approximation relies on the fact that for r ≪ R * one can safely neglect the centrifugal potential, the asymptotic kinetic energy and hyperfine splittings. On the other hand, it also ignores the higher order dispersion terms C 6 /r 6 , C 8 /r 8 , · · · , which only give relatively small corrections for the specific atom-ion system considered in this paper [10]. Hence, for R 0 r ≪ R * the single-channel wave function is given by the linear combination of the two WKB solutions (22).
Based on the arbitrariness of the reference potentials, we can choose values φ i at our convenience. For the following choice of short range phases the quantum-defect matrix takes a particularly simple form in the molecular (IS) basis where S(β) = 0, 1 denotes the total electron spin in the channel β, and a 0 = a s and a 1 = a t are the singlet and triplet s-wave scattering lengths, respectively. By using the unitary transformation (39) we find the quantum-defect matrix Y in the basis of the asymptotic channel states In the presence of an external magnetic field B, the total transformation matrix is a product of U and an additional unitary matrix Z(B) relating the bare (B = 0) and dressed (B = 0) channel states |Ψ β (B) : In this case the quantum-defect matrix reads Apart from the shifts of the internal hyperfine states of the atom and ion, the magnetic field in general also affects the motion of the ion as a charged particle. This leads the Landau quantization of the ion motion [35] and can result in several interesting scattering effects resembling the collisions in the quasi-1D confinements [36]. There are, however, two limits when the ion cyclotron motion can be neglected. First, when the characteristic size of the Landau orbit a L = /m i Ω with Ω = eB/m i and e denoting the ion charge, is much smaller than R * . Second, in the case when the ion is confined in the RF trap, which is typically much tighter than the size of the Landau states. This second case is however beyond the scope of the present analysis, and in the following we assume the former condition to apply.

Corrections to quantum-defect matrix for higher partial waves
According to the arguments presented in section 4, at distances R 0 ≪ R * one can expect the dependence of the quantum-defect matrix Y on ℓ to be negligeable. However, as we will show later with the example of a Na atom and a Ca + ion, this approximation works well only for the few lowest partial waves. Already for ℓ 4 one starts to observe some deviations of the quantumdefect parameter from the exact numerical solutions for a realistic potential. In order to improve the accuracy of the quantum-defect model we introduce a correction of the short-range phase φ for nonzero relative angular momenta ℓ. Such correction can be obtained within the semiclassical theory, since the modification of φ arises at distances where the real potential W ii (r) differs from V i (r), a region well described in the semiclassical approximation. We start from the semiclassical formula for the radial wave function f i (r) calculated in potential W ii (r) is the WKB phase, and r T is the classical turning point. Here, we neglect the contribution from the asymptotic kinetic energy, taking E = E ∞ i . Eq. (48) can be rewritten in the following way where . At sufficiently large distances, where the potential takes its asymptotic form (3), but still within the semiclassical regime r ≪ R * , we require that the phases of f i (r) and of the solutionf i (r) calculated for V i (r) be equal. At large r, ∆ i (r) → ∆ i (∞) = const, and this happens already within the semiclassical domain r ≪ R * , since the potentials W ii (r) and V i (r) have the same long-range asymptotics. The ℓ-dependent shift can be obtained by comparing the ∆ i (∞) computed for ℓ = 0 and ℓ > 0 In this formula, ∆ i | ℓ is calculated for the actual value of ℓ, while ∆ i | ℓ=0 represents a shift that can be incorporated into the definition of the short-range phase φ. By treating the centrifugal barrier in the expressions for k i (r) and k i (r) as a perturbation to the interaction potential, we find a particularly simple result containing partial-wave quantum number ℓ only in the prefactor of an integral where

Radiative charge transfer
In this section we develop a quantum mechanical description of the radiative charge exchange process A + B + → A + + B + hν with emission of a photon of frequency ν. Since the probability of charge exchange in heteronuclear atom-ion collisions is small, the quantum transition rates can be described in the distorted wave Born approximation (DWBA) [37] treating the molecular dynamics exactly and the interaction with the radiation field as a perturbation. If weak dipolar interactions are neglected, the molecular dynamics conserves the quantum numbers ℓ and m ℓ of the orbital angular momentum of the atomic fragments both in the initial and final molecular states. Radiative transitions are induced to first order by the electric dipole coupling, provided that the variation |∆ℓ| = 1.
We ignore at first nuclear spins, such that the collision only involves two coupled molecular channels. After angular integrals have been performed analytically using standard angular momentum techniques, the radiative charge-exchange decay rate reads (see, e.g., Ref. [38]): This equation is expressed in terms of reduced dipole matrix elements for free-free E ′ ℓ ′ |d(r)|Eℓ transitions, in which the atom and the ion remain unbound, and free-bound v ′ ℓ ′ |d(r)|Eℓ transitions, in which the colliding pair forms a molecular ion.
The |Eℓ scattering state is energy-normalized with incoming wave boundary conditions in the entrance potential of AB + state. The primed quantities |E ′ ℓ ′ and |v ′ ℓ ′ represent respectively energy-normalized scattering states with outgoing wave boundary conditions and unit normalized bound states of the A + B molecular complex. Energy conservation requires hν = E − E ′ + δ and hν v ′ = E − E v ′ + δ for free-free and free-bound transitions, where δ is the difference of ion and atom ionization potentials.
While this approach is physically transparent and has the advantage of providing detailed information on the products of the charge-exchange process, the generation of dipole moments as a function of photon energy and of free-bound matrix elements for all ro-vibrational levels is computationally costly. A significantly simpler approach discussed in [39] approximates the sum over all continuum and bound transitions with a simple average of a space-varying Einstein coefficient for spontaneous emission over the initial scattering wavefunction Eq. (54) is inspired by the exact quantum-mechanical closure relation over scattering and bound states in the final channel, and by the classical Frank-Condon principle relating the potential energy difference with the energy of a photon emitted at interatomic separation r, hν(r) = W i (r) − W f (r) [39]. In the next subsection we cast this equation in a particularly expressive form using a combination of semiclassical and MCQDT approaches.

Semiclassical model
The radiative charge transfer takes places at relatively short distances, where the electronic wave functions of the atom and ion start to overlap and the dipole matrix element is non negligible. At such distances one can calculate the matrix element ofĀ(r) in the semiclassical approximation Here, R min denotes the classical turning point in the potential U(r) of the entrance channel, and is the classical velocity of a particle in the entrance channel at the distance r. The MQDT function C −2 (E, ℓ) provides the proper scaling of the semiclassical wave function at short distances, with respect to its long-range asymptotic behavior, given by the energy-normalized functions |Eℓ . The total cross section for the charge exchange can be calculated from [40] σ tr = 2π k 2 hA(E).
In the regime of ultracold energies one can neglect the contribution of the energy E to v(r) at distances where the radiative charge transfer occurs. Similarly, one can also omit the contribution from the centrifugal barrier. In this approximation v(r) ≈ −2U(r)/µ, and the energy and angular momentum dependence enters only through the MQDT function C −2 (E, ℓ) Here, P tr = 2 ∞ R min drĀ(r)/v(r) is the probability of the photon emission during a single collision, which in our approximation is a constant. This way the cross section for the charge transfer event σ tr can be written as It is instructive to investigate the classical, high energy limit of the charge transfer rate, where we expect the Langevin theory to be applicable. Starting from Eq. (54) with the semiclassical approximation (56), then setting C −2 (E, ℓ) ≈ 1 at E ≫ E * , and replacing summation over ℓ by an integration over impact parameter b = (ℓ + 1 2 )/k, we obtain [16] σ class tr (E) = 2π 2µ Eq. (61) can be further simplified, first by replacing the inner integral by a constant probability P tr , and then by performing the remaining integral over b. We perform the integration over b to some maximal value of the impact parameter b max , given by the height of the centrifugal barrier E max (ℓ) = 1 4 E * ℓ 2 (ℓ + 1) 2 . Thus, the particle can penetrate the inner part of the potential only for E > E max , which imposes the upper bound on b: b max = 2R * /k. This yields for the cross section and for the charge transfer collision rate K tr = σ tr (E)v.
At large energies the charge transfer rate acquires a constant value. The latter formula agrees with predictions of the Langevin theory, assuming that all the classical trajectories that fall down on the scattering center lead to a reaction [22]. In our case the probability of a reaction in a single collision event is given by P tr , whereas in the case of homonuclear collisions, where the resonant charge transfer takes place P tr = 1 2 [2]. It is convenient to express the charge transfer rate in the quantum regime in terms of K class tr : In this form the charge transfer rate is a product of a constant classical rate, that depends on the particular atom-ion system, and a quantum factor Q(E) that is universal and depends only on the specific atom-ion combination through the characteristic energy E * and the short-range phase φ. In this way the whole energy dependence will be given by Q(E). In the limit of high energy (E → ∞) Q(E) → 1. In the limit of small energies (E → 0), only s wave contributes and Q(E) → 1 2 (1 + cot 2 φ). In Fig. 12 we show Q(E) averaged over thermal distribution for selected values of the short range phase. We observe that even at large energies E ∼ 10 4 E * , Q(E) = 1, due to the contribution of shape resonances. The detailed structure of resonances, however, is washed out by the thermal average.

Results for 40 Ca + and 23 Na
We now apply our quantum-defect model to 40 Ca + and 23 Na, a system whose ab initio potential energy curves are known [16]. This allows us to compare the predictions of our analytical approach with the full solution obtained by solving numerically the coupled-channel Schroedinger equation (1) with potentials of Ref. [16].  Table 2. Scattering channels of 40 Ca + and 23 Na for M F = 1/2 in the asymptotic representation.

Channel states
The 40 Ca ion has vanishing nuclear spin (i 1 = 0, s 1 = 1/2). We label its ground-state sublevels in an external magnetic field in increasing energy |a 1 ≡ |f = 1/2, m f = −1/2 and |b 1 ≡ |f = 1/2, m f = 1/2 . The 23 Na has nuclear spin i 2 = 3/2, hence its hyperfine angular momentum can take the values f = 1 and 2. The hyperfine structure of 23 Na is shown in Fig. 13. The figure adopts the standard |a 2 , |b 2 , · · · notation to label states in increasing energy in the magnetic field, where in the weak field limit one has the identification |a 2 ≡ |f = 1, m f = 1 , |b 2 ≡ |f = 1, m f = 0 , · · · . In the presence of an external magnetic field, only the projection M J of the total angular momentum J = f 1 + f 2 + l is conserved during a collision. However, if we ignore small anisotropic spin-spin interactions giving rise to the long-range dipole-dipole force, the states with different angular momentum l will not be coupled, and both ℓ and m ℓ will be conserved. Hence, M F will also be conserved and we can restrict our discussion to subspaces of constant M F . In most of our calculations we will consider collisions within the subblock of M F = 1/2, which contains the lowest energy channel state |a 1 a 2 . For M F = 1/2 there are four possible scattering channels, listed in Table 2 together with their threshold energies at B = 0. Table 3 lists the channel states in the molecular basis. The quantum-defect matrix Y in the molecular {IS} basis is diagonal and for the assumed  channel numbering it reads Using the frame transformation (39) one can easily find Y in the asymptotic channel representation where 1/a c = 1/a s −1/a t is the coupling parameter characterizing the strength of coupling between channels. We note that for similar triplet and singlet scattering lengths 1/a c = 0, the channels are uncoupled, and no interchannel resonances occur.
In the presence of a magnetic field, the energies of the hyperfine states of Na atom are given by the Breit-Rabi formula [41] E (Na) where g is Landé factor, and µ B is the Bohr magneton. By including in the Hamiltonian additional small couplings of the magnetic field to the nuclear spin, H(B) = H(0) + µ B B(g j j + g i i), the effective g-factor is given by g = g j − g i , and the levels acquire an additional B-dependent shift ∆E (Na) f m f (B) = g i µ B Bm f . Here, g j and g i are Landé factors for the total orbital angular momentum j = L + s and the nuclear spin i, respectively. In Ca ion, which has no hyperfine structure, the levels shift according to the standard Zeeman formula The transformation matrix Z from the bare (B = 0) to the dressed channel states (B = 0) can be easily found from the transformation matricesZ (k) , connecting bare and dressed states of ion (k = 1) and atom (k = 2). For the ion we have simplyZ (1) = 1, while for the atom wheref 2 = i 2 + 1/2 = 2, and w = y 2 + (x + Ef 2 m f2 ) 2 . Then the total transformation matrix Z for the assumed channel numbering reads

Bound states
In this section we discuss properties of bound states in a single scattering channel. Fig. 14 shows a sample bound state spectrum for a 3 Σ + molecular potential. The numerical calculations based on the ab initio Na-Ca + potentials of Ref. [16] with the quantum-defect theory with and without the ℓ-dependent correction (51) to the short-range phase are compared. For simplicity in the numerical calculation we neglect the higher-order dispersion terms in the potential, so to isolate the effect of the centrifugal barrier, in particular for high order partial waves. We note that the quantum-defect model assuming the same short-range phase φ for all partial waves starts to deviate already for For ℓ = 4 still some small discrepancy between the numerical and corrected MQDT results can be observed. Tha latter may originate from effects going beyond the WKB approximation, which was assumed to hold in the derivation of Eq. (51).

Radiative charge transfer
In the case of Na-Ca + collisions the charge transfer process occurs due to transitions from a singlet A 1 Σ channel to the absolute ground molecular term NaCa + (A 1 Σ) → Na + Ca(X 1 Σ). For its description we first use the DWBA leading to Eq. (53), which allows the final states of the charge exchange process to be individually identified. In order to calculate the free-free transition dipole elements a two coupled channel model is setup, comprising W i (r) = V (A 1 Σ) and W f (r) = V (X 1 Σ) diagonal molecular potentials with ℓ and ℓ ′ centrifugal barriers, coupled by a ξd(r) term, with ξ a parameter to be optimized. In a field-dressed approach, a vertical energy shift is introduced so to guarantee total, photons plus atoms, energy conservation. If the ξ parameter is chosen so to insure validity of first order perturbation theory, in the DWBA the required reduced matrix element is simply proportional to the transition scattering matrix element s ξ = s(E ′ ℓ ′ ← Eℓ; ξ) obtained numerically The free-bound reduced dipole elements are obtained following a different method. As a first step, the vibrational wavefunctions of the X 1 Σ potential are generated for each angular momentum ℓ ′ using a grid method. To this aim use of the step-adaptive approach of Ref. [42] is essential. In fact, due to the long-range character of the r −4 interaction, numerical convergence of near threshold levels is attained for a grid extending to large distances. To keep the number of grid points reasonable (∼ 5000), the local step size is made to vary over five orders of magnitude from the potential well to the asymptotic region using an appropriate scaling function [42]. Fig. 15 shows for instance the last rotationless vibrational wave function for a = 5R * , as well as the convergence rate of the corresponding eigenenergy for two different step sizes δr of the numerical grid near the potential minimum. Note that for this specific molecular level, convergence begins to be attained for a box of 10 5 a 0 size. The relative energy accuracy is on the order of 10 −2 for δr = 0.2a 0 , and increases by one order of magnitude for δr = 0.1a 0 . As a second step, the overlap v ′ ℓ ′ |d(r)|Eℓ with the initial scattering wavefunction of angular momentum ℓ is computed using a standard propagation code. The number of needed ℓ values is determined by numerical convergence.
In the alternative approach based on the approximate sum rule Eq. (54), the expectation value Eℓ|Ā(r)|Eℓ are extracted numerically from the elastic scattering matrix element s ξ = s(Eℓ ← Eℓ; ξ) for a single channel potential W i with angular momentum ℓ perturbed by a ξĀ(r) term. The scaling parameter ξ has to be chosen so to be in the linear regime. With this proviso, letting s 0 be the scattering matrix for the unperturbed W i potential, the needed matrix element is determined from the DWBA as The prefactor s * 0 arises from the fact that the expectation value on the right hand side has to be evaluated between scattering states with outgoing boundary conditions whereas transition matrix elements in the DWBA present incoming wave boundary conditions in the exit channel.
The charge transfer rates we obtain are shown in Fig. 16 as a function of the collision energy for a sample value of the singlet scattering length, a s = R * . The figure compares the numerical result calculated by summing contributions from all free-free and free-bound transitions Eq. (53) with the MQDT model assuming the semiclassical charge-transfer probability Eq. (64). In the former case we additionally plot separately contributions from free-free and free-bound transitions. The MQDT calculation includes the ℓ-dependent correction (51). We observe that at energies larger than 10µK, the rate of charge exchange exhibits several peaks due to the shape resonances. The MQDT model agrees well with the full numerical calculations, except in the mK regime where high order partial waves are contributing. In this range of energies the discrepancy is due to the corrections to the short-range phase for large ℓ that are beyond the applicability of the semiclassical formula (51). Finally, the approximation based on Eq. (54) is in full agreement with the numerically exact result and is not shown in the figure.

Population of vibrational states in the charge transfer process
According to Fig. 16 the main outcome of the charge transfer process are molecular ions. It is interesting to analyze the dependence of the charge exchange rate on the vibrational quantum number of the molecule, i.e. the vibrational distribution of the product molecular ions. In the full quantum treatment based on the numerical calculation of initial and final wave functions, such distribution is proportional to the transition matrix elements in Eq. (53). In the MQDT approach it can be obtained using the semiclassical approximation, in analogy with the derivation presented in Section 5.1. We start from Eqs. (54) and (56), and obtain the probability of the charge transfer process normalized per vibrational quantum number ϑ in the final state Here, ρ(E) = dϑ dE is the density of states in the exit channel and V (r) = W i (r) − W f (r) is the difference between interaction potentials of the entrance (initial) and the exit (final) channels. The distance r at which all r-dependent quantities are calculated is related to the final energy by the Frank-Condon principle: of the dissociation energies of the entrance and exit channels, respectively. The density of states ρ(E) can be calculated from the LeRoy-Bernstein formula [43]. In the the case of r −4 potential it yields In Fig. (17) we show the charge transfer rate from A 1 Σ + channel to the bound states of the X 1 Σ + channel, calculated numerically in the DWBA. The numerical points are obtained for two values of the singlet scattering length: a s = ±R * . The solid line shows prediction of Eqs. (75), with the total cross section related to the probability A by Eq. (58). We note that the full quantum result predicts oscillations of the distribution over final states, whereas the semiclassical theory leads to a smooth behavior, which can be interpreted as the distribution averaged over the quantum oscillations. The origin of these oscillations in the quantum result is due to the effects of the phase matching between the entrance and the exit channel wave functions, which is not present in the semiclassical result assuming Frank-Condon approximation.

Magnetic Feshbach resonances at zero energy
In this section we analyze magnetic Feshbach resonances in the limit of zero energy when only s-wave scattering is present. We neglect the effects of the magnetic field on the translational ion motion. For the Na-Ca + system of interest, the radius of the Landau orbit a L becomes equal to R * at B L = 1086G. Hence, the presented analysis is valid as long as B ≪ B L . Figure 18 shows the variation of the s-wave scattering length versus magnetic field for two sample values of the singlet and triplet scattering lengths: a s = R * = 2081a 0 and a t = −R * = −2081a 0 , respectively. The resonances can be assigned to the particular bound state crossing the threshold of the α = |a 1 a 2 channel, as one can verify by inspection of the bound state energies in the bottom panel. The resonances are labeled by capital letters from A to E. The middle panel shows a close-up of the bound-state spectrum just below threshold. Because of the long-range character of the polarization potential, the last s-wave bound state is always located relatively close to the threshold, with the binding energy E b ≤ 106E * = 3.02 MHz. As one can observe in the figure the last bound state in the open channel can be strongly coupled to the other bound states crossing the threshold, giving rise to relatively strong avoided crossing, as in the case of resonance E. In contrast, very close to threshold the energy of the bound state bends and follows the universal behavior E = − 2 (2µa 2 ).
An approximate but highly accurate description of the Feshbach resonances can be developed using a two channel description based on CI (configuration interaction) model of a Feshbach resonance, where a single closed channel represents the effects of all the closed channels contributing to a resonance [44]. In the case of a single open and a single closed channel, the application of MQDT is straightforward, and after some simple algebra (see for instance [45]), one obtain the following expression for the phase shift ξ(E, l) in the open channel [46] ξ . (77) The first term ξ bg (E, l) is the background phase shift, describing the scattering from the open channel only and incorporating effectively the influence of the closed channels for magnetic fields far from the resonance. The second term describes the resonant contribution resulting from a bound state in the closed channel with energy E n , crossing the threshold of the open channel. An energy-dependent width of the resonance is given by a constant width Γ multiplied by the MQDT function C −2 (E, ℓ), which accounts for a proper threshold behavior as k → 0. In the two-channel model one assumes that the energy of a bound state varies approximately linearly with the magnetic field B: where B n is the magnetic field at which the bound state crosses the threshold of the open channel, and δµ is the difference of magnetic moments between the open and closed channel. The second MQDT function tan λ(E) in the denominator of ξ res , describes the shift of the resonance position from the bare value B n , which is due to the coupling between the open and closed channel. In the following we will focus on s-wave Feshbach resonances, considering the zero-energy limit. According to the results of Sections 3.2 and 3.4, for ℓ = 0 the phase shift and MQDT functions exhibit the following threshold behavior where s = a/R * . The connection to the standard theory of magnetic Feshbach resonances can be done by introducing the width of magnetic Feshbach resonance ∆: The parameter ζ describes the range of magnetic fields expressed as a fraction of the resonance width ∆ over which the resonance exhibit the universal properties and the occupation of the closed channel remains small. Sufficiently close to the resonance, the binding energy is given by the universal formula E = − 2 (2µa 2 ), which leads to [47] For entrance channel dominated resonance ζ ≫ 1 and Z(B) remain small for detuning of the order of ∆. In the opposite case ζ ≪ 1, the Feshbach resonance is called closed-channel dominated, and the universal regime where the energy-independent formula (84) is applicable is very narrow. We have fitted the universal formula (84) to Feshbach resonances presented in Fig. 18. The results are summarized in Table. 4. The difference of the magnetic moments has been determined from the bound state spectrum (lower panel of Fig. 18) by extracting the linear slope for the molecular states giving rise to the resonances. Finally, the parameter ζ has been calculated from Eq. (87). We note that only the two first Feshbach resonances, occurring at relatively small magnetic fields, are non universal, while the remaining ones are relatively broad and entrancechannel dominated.

Magnetic Feshbach resonances at finite energies
We analyze here the dependence of elastic and inelastic rates on the magnetic field at finite collision energy, where several partial waves play in general an important role in the collision physics. To this aim we use both our MQDT model and numerical calculations based on the close-coupled (CC) Schroedinger equation. In order to calculate the charge exchange rates from the CC solutions   Fig. 18: resonance position B 0 , resonance width ∆, background scattering length a bg , difference of magnetic moments between the open and closed channel δµ, and parameter ζ characterizing the fraction of the closed channel (see Eq. (87)). The labels in the first column enumerate the resonances in accordance with markings in Fig. 18. we have generalized the approximate closure relation Eq. (54) to the multichannel case, including hyperfine effects and the external magnetic field. The charge-exchange process in the MQDT model has been described using the semiclassical approximation Eq. (60) where the MQDT function C −2 (E) is replaced by the amplitude of the wave function in the singlet A 1 Σ channel. Details are provided in Appendix Appendix B. Figs. 19, 20 and 21 show the elastic and charge transfer rates at three collision energies: E = 1nK, E = 1µK and E = 1mK for some typical values of the singlet and triplet scattering length, a s = R * = 2081a 0 and a t = −R * . Contributions from the lowest 2, 9 and 21 partial waves are included respectively at these three energies. At 1nK resonances appear mostly in the s-wave and are relatively broad. In addition one can observe few narrow resonances in the p-wave channel occurring at B = 27, 84.7 and 182G. In this range of temperatures the MQDT model is extremely accurate and agrees perfectly with the full numerical CC calculations.
At higher energy E = 1µK the agreement of the analytic model with the numerical solution is still very good. At this energy the resonance peaks arise from the 4 lowest partial waves. The arrow in the bottom panel of Fig. 20 indicates the ℓ = 3 Feshbach resonance. It is the only resonance at E = 1µK whose position is not well predicted by the analytical MQDT model. We note that the charge transfer rates exhibit more resonance peaks than the elastic rates. This general behavior can be qualitatively understood by analyzing the number of partial waves contributing to the elastic and inelastic processes. The elastic rates are dominated by the reflection from the long-range polarization potential and their contribution to the cross section decay as (2ℓ + 1) sin 2 ξ ℓ (k) ∼ 1/ℓ 5 at large ℓ (see Eq. (33)). In contrast the charge exchange process must involve the tunneling through the centrifugal barrier, and its probability decays exponentially with ℓ. Hence the number of partial waves contributing to the inelastic process is much smaller than for elastic scattering, and the inelastic rates are more sensitive to scattering resonances. On the other hand, the narrow resonances from high-order partial waves that appear in the elastic rates are less pronounced due to the strong background arising from reflection on the long-range potential.
Finally, at the highest energy considered E = 1mK the resonances are narrower and have smaller amplitude. This effect arises due to the large number of partial waves contributing to the scattering, which have tendency to wash out the resonance structures. We observe that at 1mK the MQDT model basically follows the magnetic field dependence of the exact numerical rates, but it predicts accurately only the resonances associated to the lowest partial waves. We note that inclusion of the thermal averaging washes out the resonance structure for the elastic rates, while the charge-transfer rates still exhibit some resonance peaks. This is again due to the vastly different numbers of partial waves contributing to the elastic and inelastic collision processes.

Conclusions
Summarizing, we have developed a quantum-defect model for ultracold atom-ion collisions. The model was applied to the reference system composed by a 40 Ca + ion and a 23 Na atom, and its predictions were thoroughly verified by comparison with numerical close-coupled calculations using ab initio potential energy curves. Our model is based on the multichannel quantum-defect formalism, where the quantum-defect parameters are defined in terms of the analytic solutions for r −4 polarization potentials. Use of a frame-transformation allows us to reduce the number of  short-range parameters to essentially singlet and triplet scattering lengths only. Since for atomion systems of experimental interest the values of the singlet a s and triplet a t scattering lengths are not yet known, in our calculations we have assumed typical scattering lengths of the order of the characteristic length R * . Once a s and a t will be measured our model could be readily applied to obtain all the basic collisional properties in the ultracold domain, including accurate positions of Feshbach and shape resonances. Application of our theory to other atom-ion systems is straightforward, amounting to a simple change of the scattering lengths and to the use of new characteristic parameters R * and E * , which are determined by the atomic polarizability and the reduced mass. In our studies we considered only two-body collisions in free space ignoring possible effects of the trapping potential, which is present in realistic systems. This can be of particular importance for the ions, that acquire in the presence of a time-dependent radio-frequency potential a small amplitude high frequency motion known as micromotion. This effect is important for current experiments, leading on one side to a significant loss of atoms, on the other potentially preventing sympathetic cooling of the ions to the ground-state of the trapping potential [48]. Moreover, the presence of a tight ion trap with characteristic size smaller than R * modifies the long-range asymptotics of the atom-ion wave functions and results in principle in the coupling of relative and center-of-mass motion. In addition, in the presence of a magnetic field the charged ions describe cyclotron orbits, an additional confinement effect which can lead to the appearance of scattering resonances [36] and may affect the actual positions of Feshbach resonances. We are currently investigating these issues.
To find values of the coefficients c n it is sufficient to set h + m = 1 and h − m = 1 for some, sufficiently large m and calculate h + n and h − n up to n = 1 using (A. 8 In numerical calculations it is more convenient to find ν from equation [32] cos πν = 1 − ∆(1 − cos π with γ n = q/(4n 2 −a). Typically determinant ∆ converges very fast and to calculate ∆ it is enough to take relatively small matrices.
Appendix B. Multichannel calculations of the radiative charge transfer Equation (54) based on the approximated closure relation can be easily generalized it to the multichannel case including hyperfine effects and an external magnetic field. In this case the initial state for atoms incoming in the dressed channel state α of Eq. (46) will be labeled |Eℓα . The reduced dipole moment becomes a matrix, with only nonvanishing diagonal elements if the total electron and nuclear spin {SI} representation is used.
In the presence of a magnetic field close-coupled scattering equations are numerically solved with a matrix perturbation ξĀ(r). The emission rate for atoms incoming at collision energy E is expressed in the DWBA in terms of the open-open scattering matrix elements s ξβα = s(Eℓβ ← Eℓα; ξ) with and without field perturbation As in the single-channel case, right multiplication by the s * 0 matrix enforces the correct boundary conditions, and symmetry of the scattering matrix (s ξαβ = s ξβα ) resulting from time-reversal invariance has been used.
Calculation of the charge-exchange rates in the MQDT approach can be done using the semiclassical formula Eq. (60), with the function C −1 (E, ℓ) replaced by the amplitude A Here, H ij (kr) → δ ij ie −i(k i r−ℓ i π/2) / √ k i and H where χ ββ ′ = δ ββ ′ φ S(β) is the diagonal matrix containing short-range phases of the singlet φ 0 and the triplet φ 1 potentials, with S(β) denoting the total electron spin in the channel β.