Mathematical Modeling of Resonant Processes in Confined Geometry of Atomic and Atom-Ion Traps

We discuss computational aspects of the developed mathematical models for resonant processes in confined geometry of atomic and atom-ion traps. The main attention is paid to formulation in the nondirect product discrete-variable representation (npDVR) of the multichannel scattering problem with nonseparable angular part in confining traps as the boundary-value problem. Computational efficiency of this approach is demonstrated in application to atomic and atom-ion confinement-induced resonances we predicted recently.


Introduction
During the last three decades one can observe an impressive development of the physics of ultracold atoms [1] and cold ions [2]. Different aspects of this investigations attract big interest from the side of theoreticians and experimentalists. Thus, experimentalists have got a chance to work here with deterministically prepared quantum systems [3] with precise control of interparticle interaction, particle states and particle number [2,3]. It opens a possibility for quantum simulation with fully controlled few-body systems in ultracold quantum gases [4][5][6]. The control of interparticle interaction is performed in quantum gases with the help of magnetic Feshbach resonances [1] which, in confined geometry of atomic traps, transform into so-cold confinement-induced resonances (CIRs) [7][8][9][10][11][12][13].
The mathematical modeling of resonant processes in the confined geometry of atomic traps is a special computational problem due to the strong nonseparability of the angular part of the system wave-function, what is the consequence of the strong particle-trap interaction. Therefore, the conventional theory is no longer valid here and the development of the low-dimensional theory, including the influence of the confinement, is needed. In our works we have developed computational methods [9,10,14,15] for pair collisions in tight atomic waveguides and have found several novel effects in its applications: the CIRs in multimode regimes including effects of transverse excitations and deexcitations [10], the so-called dual CIR yielding a complete suppression of quantum scattering [9], and resonant molecule formation with transferring energy relies to center-of-mass excitation while forming molecules [16]. The last effect was recently confirmed in the Heidelberg experiment [3]. Our calculations have also been used for planning and interpretation of the Innsbruck experiment on e-mail: melezhik@theor.jinr.ru

Resonant Scattering in Atomic and Atom-Ion Traps
The mathematical modeling of resonant scattering in confined geometry of atomic and atom-ion traps requires the development of efficient computational methods for integration of the system of 3D Schrödinger-like equations strongly coupled by the matrix of effective interatomic potentialV(r) (Î is unit matrix here) and the derivation of the wave-function satisfying the scattering asymptotic at |z| → +∞ for the fixed collision energy E, where r = (x, y, z) = (ρ cos φ, ρ sin φ, z) = (r sin θ cos φ, r sin θ sin φ, r cos θ) is the relative variable between the colliding atoms [9][10][11][12][13] or the atom variable relative the heavy ion [18]. Here, f e (E) is the desired scattering amplitude, Φ 0 (x, y) is the wave function of the groundstate of the 2D harmonic oscillator 1 2 µ(ω 2 x x 2 + ω 2 y y 2 ) describing the interaction of the atoms with the confining trap and k 0 = 2µ(E − ω ⊥ )/ = 2µE / → 0 is the relative momentum of two atoms in the open channel e, where ω ⊥ = (ω x + ω y )/2 and µ is the reduced mass of two colliding atoms [9][10][11][12][13] or the mass of atom colliding with heavy ion [18].
In the case of the atom-atom interaction a four-channel potentialV(r) describing three lowest magnetic Feshbach resonances in Cs [19] (with c = 1, 2, 3) was used in works [11,12]. Some of these resonances have been used in experimental investigation of the Cs confinement-induces resonances in anisotropic waveguide-like traps [8]. In the case of atom-ion confined scattering the problem (1, 2) was reduced to one channel scattering (α = e) with scalar atom-ion potentialV(r) → V(r) ∼ C 4 /r 4 [18].
Note, that in the presence of an anisotropic harmonic trap (ω x ω y ) the problem (1, 2) becomes non-separable in the 3D space {ρ, z, φ}. To resolve the problem in the case of atom-atom confined scattering [11,12] we have to integrate the system of four coupled 3D Schrödinger-like equations (1). In an isotropic trap (ω x = ω y ) the problem admits separation of the angle φ and reduction to 2D case. In free space (ω x = ω y = 0) the angular part separates and the problem reduces to four coupled radial equations considered in [19]. The choice of the potentialV(r) in the tensorial form (α = {e, c}) with the diagonal terms V αα (r) as square-wells enabled the derivation [19] of an analytic model for describing Feshbach resonances in free space (ω x = ω y = 0). However, to develop an efficient computational scheme for the integration of the problem (1, 2) in a confining trap (ω x 0, ω y 0) is a challenging task due to sharp jumps at the edge of the interaction in the diagonal terms V αα (r) chosen in the form of square-well potential [11,12]. In the case of atom-ion confined scattering the tensorial structure of the potentialV(r) is absent but the need arises to get accurate solution for the long-range scalar atom-ion interaction V(r) ∼ C 4 /r 4 [18].

Method
The key feature of the system of equations (1) is the strong nonseparability of the angular part of the desired wave-function |ψ(r) due to the strong coupling of the partial waves defined by the atomtrap interactionŴ(r) = 1 2 µ(ω 2 x x 2 + ω 2 y y 2 ). To solve the difficulty, the approximation of the angular part is done within the nondirect product discrete-variable representation (npDVR) which was suggested [20][21][22][23] and developed in our works [9-13, 16, 22-31]. It was shown that 2D npDVR is much more efficient computationally as compared to the conventional partial wave analysis due to its fast convergence and flexibility: there is no need for laborious calculation of the matrix elements under the change of the form of the interactions because any local interaction is a diagonal one in the npDVR.
In npDVR we have formulated the multichannel scattering problem with nonseparable angular part as a boundary-value problem and successfully applied this scheme to scattering in confined geometry of atomic and hybrid atom-ion traps [10][11][12][13]18]. Here we discuss computational aspects of this approach and recent results obtained in quantitative analysis of the atom-ion CIRs.

Nondirect Product Discrete-Variable Representation
Different kinds of the 1D DVR or Lagrange-meshes are widely used in quantum mechanics computations [32][33][34][35] due to the simplicity (DVR leads to diagonal form for the interaction and compact form for the kinetic energy operator) and efficiency (fast convergence and stability) of this approach. In the work [20] we have theoretically estimated and confirmed by numerical experiment that the error very quickly decreases ∼ 1/(N + 1)! with an increase in the number of DVR basis functions N. It ensures fast convergence of the 1D DVR. However, an extension of this representation to the 2D case (two angles θ and φ of the unit sphere) is a nontrivial problem. Actually, the simple idea to construct the 2D DVR as a direct product of two 1D DVRs leads to essential complication of the matrix of the angular part of the kinetic energy operator [33,36]. As a result, the advantages of the 1D DVR, its simplicity and efficiency, are lost [33]. An alternative way to construct 2D DVR on an unit sphere is to use the spherical harmonics defined on a two-dimensional grid. However, in this case it becomes not possible to satisfy the orthogonality conditions for all the elements of the fixed set of this basis on the chosen grid [36,37]. To overcome this difficulty we have suggested [22,23] to use the basis of the orthogonalized combinations of the spherical harmonics on the two-dimensional grid over θ and φ variables. It turned out that this idea was very efficient for the time-dependent, as well as, stationary Schrödinger equations with two [10][11][12][13]18], three [22][23][24][25][26][27][28] and four [9,16,29,30] nonseparable spatial variables. Our idea how to construct the 2D npDVR [22,23] was also successfully extended for computing vibrational levels of four-atom molecules [36] where it has got the title non-direct product DVR (npDVR). An alternative variant of 2D npDVR with non-direct product angular grid coinciding on the unit sphere with the nodes of the Lebedev quadratures was suggested in [37].
By applying the 2D npDVR for Eq. (1) we get the system of 4N Schrödinger-like coupled equations with respect to the unknown vector u(r) = { λ j u α j (r)} (here j = 1, . . . , N; α = e, c) where Here x j = r sin θ j cos φ j and y j = r sin θ j sin φ j and the elements u α j (r) of the vector u(r) coincide with the values rψ α (r, Ω j ) of desired wave function on N angular grid points Ω j = {θ j , φ j }. The matrix Y jν and its inverse (Y −1 ) ν j are defined on the grid Ω j as described in [9,10,15,24] (here one can also find the procedure of construction of the 2D angular grid Ω j and the definition of the weights λ j ).
Thus, by using the 2D npDVR we transform the initial scattering problem (1, 2) to the system of the Schrödinger-like coupled equations (3) with the following asymptotic u e j (r) = r e ik 0 z j + f e e ik 0 |z j | Φ 0 (x j , y j ) , u c j (r) = 0 at |z j | = |r cos θ j | → +∞ and u α j (r) = 0, at r = 0 , which follows from the finiteness of the desired wave function ψ α (r) (α = e, c) at r → 0. = r m−1 cos θ j (let say for j θ ≤ j m and for j θ ≥ N θ − j m + 1, here j m is a constant) fall into the asymptotic region |z m j | = |r m cos θ j | → +∞ where the wave function u e j (r) satisfies the asymptotic formulae (7) (see Fig. 1).

Boundary conditions: reduction scattering problem to boundary-value problem
where x m j = r m sin θ j cos φ j and y m j = r m sin θ j sin φ j . Thus, the initial scattering problem (1, 2) is reduced to the boundary-value problem for the system of ordinary differential equations of second order (3) with the boundary conditions (8) and (9) which do not contain the unknown scattering amplitude f e . After integration of the problem (3), (8) and (9) one can find the scattering amplitude f e by mapping the calculated wave function ψ α (r, Ω j ) at the points |z m j | = |r m cos θ j | → +∞ with the asymptotic boundary condition (7). We solve the system of Eqs. (3) on a quasi-uniform radial grid [22,27] r n = r m e γζ n − 1 e γ − 1 , n = 1, 2, . . . , N r (10) of N r grid points {r n } defined by mapping r n ∈ [0, r m ](r m → +∞) onto the uniform grid ζ n ∈ [0, 1] with the equidistant distribution ζ n − ζ n−1 = 1/N r . One can achieve a suitable distribution of the grid points for a specific interatomic and confining potential by varying N r and the parameter γ > 0.
To integrate the boundary-value problem (3), (8) and (9) one can apply efficient algorithms. By using high order (six-order) finite-difference approximation for the radial derivatives on quasi-uniform grid we arrive at a system of algebraic equations with block-band structure which can be solved by the LU decomposition [38] or the sweep (also known as the Thomas algorithm [39]) method [40], which are very efficient for speed and memory.

Results
The high computational efficiency of the mathematical model (3,8) and (9) was demonstrated in application to the resonant scattering in confined geometry of waveguide-like atomic traps [10][11][12]. In Fig. 2 we present calculated atomic probability density distribution |ψ e (x, z)| 2 in the vicinity of d-wave Cs Feshbach resonance (occurs at B 0 = 47.8 G [19]) for B = 47.9 G in the trap with ω ⊥ = ω x = ω y = 14 kHz [8]. This resonance is developed as a peak at the region r a of the action of the square-well potential between the atoms in the plot of probability density distribution with a rapid oscillating behaviour over r (which is not clear in the figure), due to the constitution of quasi-molecular states. In the asymptotic region |z| → +∞ we observe slow oscillations of the channel wave function ψ e (x, z) (2) with the period of oscillations 2π/k 0 → ∞ defined by the very small colliding energy k 0 = 2µE / → 0. Fig. 3 illustrates the convergence of the method with respect to the number of grid points N over the angular variable on a typical example of calculating the T coefficient T = |1 + f e | 2 near s-wave Cs Feshbach resonance (which occurs at B 2 = −11.16 G [19]) for B = 14.14 G. It demonstrates rather fast convergence of the 2D npDVR over angular grid points N despite the strong coupling over angular variable provided by the atom-trap interaction potential 1/2µω 2 ⊥ r 2 sin 2 θ in a wide area of variation of the interatomic distance r. It is shown that for getting absolute accuracy on the level 10 −3 it is enough to keep about N = 10 of equations in the system (3). Further computations with increasing N demonstrate monotonic convergence to more accurate T coefficient. To reach the level of accuracy 10 −4 it needs to increase N up to 30. With this method the atom-ion CIRs were predicted and quantitatively investigated [18]. In Fig. 4 we present calculated transmission coefficients T (R * /a ⊥ , a ⊥ /a 3D ) for the atom-ion scattering in confined geometry of hybrid atom-ion trap ω ⊥ = ω x = ω y as a function of ratios of the range of atom-ion interaction R * to the width of atomic trap a ⊥ = /(µω ⊥ ) and a ⊥ to the atom-ion scattering length a 3D in free space. The minima in the coefficient T indicate the positions of the atom-ion CIRs. Before this investigation the CIR were known only for confined atomic scattering at the point a ⊥ /a 3D = 1.46 for R * /a ⊥ → 0 (see arrow in Fig. 4) [7][8][9][10][11][12][13]. Our investigation [18] has predicted conditions for CIR appearance in the confined atom-ion scattering. It was shown that due to longer tail of the atom-ion interaction with respect to atom-atom interaction, the limitation R * /a ⊥ → 0 characterizing the atomic 6 EPJ Web of Conferences 173, 01008 (2018) https://doi.org/10.1051/epjconf/201817301008

Mathematical Modeling and Computational Physics 2017
CIR is overcame in hybrid atom-ion confined systems. Due to the importance of the CIRs for controlling the interparticle interactions further investigations are in demand. Specifically, the effect of the ion-micromotion induced heating is most actual unsolved problem. It asks for additional development of efficient computational methods [2]. Our computational schemes are look very promising in this respect. The transmission coefficients T (R * /a ⊥ , a ⊥ /a 3D ) calculated for the confined pair Li-Yb + collisions for different values R * /a ⊥ and a ⊥ /a 3D .

Conclusion
In this lecture we have considered computational aspects of the mathematical models we developed recently for ultracold two-body resonant processes in atomic and hybrid atom-ion traps. The main attention is paid to formulation in the npDVR of the multichannel scattering problem with nonseparable angular part in the confining traps as a boundary-value problem. The computational efficiency of this approach is demonstrated in application to atomic and atom-ion confinement-induced resonances we predicted recently.
The computational efficiency of our mathematical models achieved in quantitative analysis of different resonant processes and the supposed transparent procedure for its parallelization make the developed techniques very promising in application to other challenging problems of low-dimensional few-body physics. Note here the computationally demanding, and actual in physics of quantum gases and ions, problem of the impact of the ion-micromotion induced heating on the surrounding atomic gas in hybrid atom-ion ultracold systems [2,18].