Asymptotic model for shape resonance control of diatomics by intense non-resonant light

To derive the asymptotic approximation, we consider the Schrödinger equation with V_g(R) replaced by its leading order asymptotic term, ${{V}_{g}}(R)\approx -{{C}_{6}}/{{R}^{6}}$ describing the van der Waals interaction. For the interaction with the non-resonant field, we also account only for the leading order term which scales as $1/{{R}^{3}}$. In addition, the R-independent term in ${{\alpha }_{\bot }}(R)$, lowers the dissociation limit. It reduces to ${{E}_{0}}=-\frac{4\pi }{c}{{\alpha }_{0}}I$, where ${{\alpha }_{0}}$ denotes the atomic polarizability. Taking advantage of m being conserved and considering only m = 0, the asymptotic 2D-Schrödinger equation reads

$$\begin{eqnarray}&&\begin{array}{l} \left[ -\frac{{{\hbar }^{2}}}{2\mu }\frac{{{{\rm d}}^{2}}}{{\rm d}{{R}^{2}}}-\frac{{{C}_{6}}}{{{R}^{6}}}+\frac{{{\hbar }^{2}}}{2\mu }\frac{{{{\bf L}}^{2}}}{{{R}^{2}}} \right]\psi (R,\theta ,\varphi ) \\ \qquad -\frac{2\pi I}{c}\left[ 2{{\alpha }_{0}}+\frac{2\alpha _{0}^{2}}{{{R}^{3}}}(3{{{\rm cos} }^{2}}\theta -1) \right]\psi (R,\theta ,\varphi )=E\psi (R,\theta ,\varphi ). \\ \end{array}\end{eqnarray} \tag{ 3 }$$

If the atoms are not identical, $\alpha _{0}^{2}$ needs to be replaced by ${{\alpha }_{1}}{{\alpha }_{2}}$, and $2{{\alpha }_{0}}$ by ${{\alpha }_{1}}+{{\alpha }_{2}}$, see equation (2). Equation (3) can be rescaled by introducing a dimensionless reduced length x, a reduced energy $\mathcal{E}$ (defined with respect to the field shifted dissociation limit E₀) and a reduced laser field intensity $\mathcal{I}$

$$\begin{eqnarray*}\begin{array}{rcl} R & = & \sigma x, \\ E-{{E}_{0}} & = & \epsilon \mathcal{E}, \\ I & = & \beta \mathcal{I}. \\ \end{array}\end{eqnarray*}$$

The unit conversion factors for length σ, energy and laser intensity β contain the information specific to the free atom pair:

$$\begin{eqnarray}&&\sigma ={{\left( \frac{2\mu {{C}_{6}}}{{{\hbar }^{2}}} \right)}^{1/4}},\end{eqnarray} \tag{ 4a }$$

$$\begin{eqnarray}&&\epsilon =\frac{{{\hbar }^{2}}}{2\mu {{\sigma }^{2}}},\end{eqnarray} \tag{ 4b }$$

$$\begin{eqnarray}&&\beta =\frac{c}{12\pi }\frac{{{\hbar }^{3/2}}C_{6}^{1/4}}{\alpha _{0}^{2}{{(2\mu )}^{3/4}}}=\frac{c{{\sigma }^{3}}\epsilon }{12\pi \alpha _{0}^{2}}.\end{eqnarray} \tag{ 4c }$$

The unit conversion factor for time, $\hbar /\epsilon$, is obtained from that of energy. The unit conversion factor for intensity, β, is proportional to $C_{6}^{1/4}\alpha _{0}^{-2}$ (i.e., to $\alpha _{0}^{-3/2}$, according to the approximate London formula for the van der Waals coefficient C₆) and to ${{\mu }^{-3/4}}$, such that larger polarizability and larger reduced mass require less intensity I for achieving the same value of the reduced intensity $\mathcal{I}$.

Employing atomic units, that is Bohr radii a₀ for σ, Hartree for and $a_{0}^{3}$ for the atomic polarizability, and expressing the laser intensity I in GW cm⁻², the reduced intensity is given by

$$\begin{eqnarray}&&\mathcal{I}=4.274177\;{{10}^{-8}}\;\frac{\alpha _{0}^{2}I}{\epsilon {{\sigma }^{3}}},\end{eqnarray} \tag{ 5 }$$

whereas the shift of the dissociation limit in reduced units is equal to

$$\begin{eqnarray}&&{{\mathcal{E}}_{0}}=\frac{{{E}_{0}}}{\epsilon }=-\frac{4\pi {{\alpha }_{0}}I}{c\epsilon }=-1.424725\;{{10}^{-8}}\;\frac{{{\alpha }_{0}}I}{\epsilon }.\end{eqnarray} \tag{ 6 }$$

When the reduced energy is expressed in μK, the numerical factors are equal to 13496.717 for the reduced intensity and −4498.93 for the reduced threshold shift, respectively. The asymptotic Schrödinger equation in reduced units is given by

$$\begin{eqnarray}&&\left[ -\frac{{{{\rm d}}^{2}}}{{\rm d}{{x}^{2}}}-\frac{1}{{{x}^{6}}}+\frac{{{{\bf L}}^{2}}}{{{x}^{2}}}-\mathcal{I}\frac{{{{\rm cos} }^{2}}\theta -1/3}{{{x}^{3}}}-\mathcal{E} \right]f(x,\theta ,\varphi )=0.\end{eqnarray} \tag{ 7 }$$

This universal asymptotic Schrödinger equation is valid at sufficiently large distances where the potential is dominated by the $1/{{x}^{6}}$ term, i.e., for $x\gt {{x}_{{\rm asym}}}={{({{C}_{8}}/{{C}_{6}})}^{1/2}}$.

The asymptotic model in reduced units predicts that a field-free shape resonance is solely determined by equation (7), i.e., by its rotational quantum number ℓ, and by boundary conditions at a distance, ${{x}_{0\ell }}\gt {{x}_{{\rm asym}}}$. The value of ${{x}_{0\ell }}$ is related to the value in reduced units of the s-wave scattering length of the molecule. In the presence of a non-resonant field, the resonance energy (in reduced units) depends, apart from the field-free scattering length, on both the reduced laser field intensity $\mathcal{I}$ and the field-free rotational quantum number ℓ.

To account for the presence of the interaction with the non-resonant light in equation (7), we introduce below a modification of the nodal line technique to solve the asymptotic Schödinger equation for $x\gt {{x}_{{\rm asym}}}$. To this end, we first expand the wave function in spherical harmonics, ${{Y}_{\ell ,m}}$ (θ, φ)(setting m = 0)

$$\begin{eqnarray}&&f(x,\theta ,\varphi )=\mathop{\sum }\limits_{\ell }{{y}_{\ell }}(x){{Y}_{\ell ,m=0}}(\theta ,\varphi ),\end{eqnarray} \tag{ 8 }$$

introducing the radial functions ${{y}_{\ell }}(x)$ for the different coupled channels ℓ. Equation (7) is then replaced by a system of coupled equations which can be written in vectorial form

$$\begin{eqnarray}&&\frac{{{{\rm d}}^{2}}}{{\rm d}{{x}^{2}}}{\bf y}(x)+({\bf M}+\mathcal{E}\;{\mathbb{1}}){\bf y}(x)=0,\end{eqnarray} \tag{ 9 }$$

where the vector ${\bf y}(x)$ is the set of functions ${{y}_{\ell }}(x)$, ${\mathbb{1}}$ denotes the identity and ${\bf M}$ is the matrix of the operator $\frac{1}{{{x}^{6}}}-\frac{{{L}^{2}}}{{{x}^{2}}}+\mathcal{I}\frac{{{{\rm cos} }^{2}}\theta -1/3}{{{x}^{3}}}$ represented in the basis of spherical harmonics with ℓ-values of the same parity. We restrict our considerations here to m = 0 and even ℓ values varying from 0 to various ℓ_max so that the model consists of $n={{\ell }_{{\rm max} }}/2+1$ channels $\ell =0,2,\ldots ,2(n-1)$. We denote by ${{{\bf y}}^{j}}(x)$ a particular solution of the asymptotic Schrödinger equation in the coupled channel model

$$\begin{eqnarray}&&{{{\bf y}}^{j}}(x)=\mathop{\sum }\limits_{\ell \,{\rm even}=0}^{{{\ell }_{{\rm max} }}}y_{\ell }^{j}(x){{Y}_{\ell ,m=0}}({\rm cos} \theta ),\end{eqnarray} \tag{ 10 }$$

where $y_{\ell }^{j}(x)$ is the radial component of the jth solution in the ℓth channel.

The nodal line technique is a method to solve the Schrödinger equation in the asymptotic region, which is applicable to the near-threshold energy range. The disregarded inner region is dominated by the potential V_g(R), cf equation (1), and gives rise to a nodal structure of the wavefunction which is almost independent of its energy and the other interaction terms in the Hamiltonian (1). Therefore the dependence of the node positions on energy $\mathcal{E}$, on the centrifugal term $\propto \;\ell (\ell +1)$, or on any other small term in the Hamiltonian, such as the coupling with the non-resonant light $\propto \;\mathcal{I}$, can be estimated in the WKB approximation by a Taylor expansion of the local de Broglie wavelength. To first order, it introduces small contributions to the node positions which vary linearly with $\mathcal{E}$, $\ell (\ell +1)$, and $\mathcal{I}$, respectively.

When using the nodal line technique, the solution of the coupled equations (9) is carried out only for the zone where the asymptotic Hamiltonian is valid. At large distance, physical boundary conditions are imposed, depending on the sign of $\mathcal{E}$. For $\mathcal{E}\lt 0$, the radial wave functions exponentially decay in all channels, quantifying the energy of bound levels. For $\mathcal{E}\geqslant 0$, regular and irregular Bessel functions characterize the asymptotic behavior. At small distance, on the frontier of the inner zone, we require the radial part of the physical wave function in each channel $y_{\ell }^{j}$ to vanish at a position ${{x}_{0\ell }}$ that is located on a ℓ-dependent straight line in the $(\mathcal{E},x)$ plane, the so-called nodal line [9, 34].

Without non-resonant field, the following positions were used:

$$\begin{eqnarray}&&{{x}_{0\ell }}={{x}_{00}}+A\mathcal{E}+B\ell (\ell +1),\end{eqnarray} \tag{ 11 }$$

where the parameters x₀₀, A and B are characteristic of the chosen atom pair. In particular, x₀₀ corresponds to the position of a node of the threshold s-wave wave function and is related to the s-wave scattering length [7]. A takes the variation of the node position with energy in the wave function with $\ell =0$ into account. B describes the shift in the node of the threshold wave functions induced by the centrifugal term for the various partial waves, $\ell \gt 0$. The parameters x₀₀, A and B are adjusted, if possible, to reproduce experimental data, such as the positions of bound levels or resonances close to threshold, and the s-wave scattering length. They can also be determined by studying the nodal structure of wave functions calculated in a single-channel model using a molecular potential, when available.

In the absence of either reliable potentials or experimental data, there is a rough, but universal estimate of the parameters A and B given by very simple analytical formulas which depend only on the s-wave scattering length, ${{A}^{G}}=-{{({{x}_{00}})}^{7}}/8$ and ${{B}^{G}}={{({{x}_{00}})}^{5}}/4$ [13]. These laws are deduced from the universal model of [4] which consists in a $-1/{{x}^{6}}$ potential limited by an infinite repulsive wall at a distance ${{x}_{{{0}_{G}}}}\to 0$. The WKB approximation is used to evaluate, in the vicinity of the threshold and for a not too high value of ℓ, the shift of the node located at x₀₀ that arises from the contribution of the kinetic $A\mathcal{E}$ and centrifugal $B\ell (\ell +1)$ energies in the range ${{x}_{{{0}_{G}}}}\leqslant x\leqslant {{x}_{00}}$ [34]. Although the model becomes less realistic as x₀₀ decreases, the corresponding A^G and B^G values are comparable to the values ajusted to experimental data [9, 10].

In the presence of a laser field, as it will be shown below (see section 4.1), an intensity dependent term has to be added to the nodal lines

$$\begin{eqnarray}&&{{x}_{0\ell }}={{x}_{00}}+A\mathcal{E}+B\ell (\ell +1)+C\mathcal{I}.\end{eqnarray} \tag{ 12 }$$

The new term, $C\mathcal{I}$, i.e., lowest order in $\mathcal{I}$, accounts for the contribution of the interaction with the non-resonant field at short range. With this modification it is possible to obtain the bound levels, the resonance profiles of the shape resonances as well as the scattering length of the field dressed molecule for any intensity. The $\mathcal{I}$-dependent term can be evaluated in the same way as A^G and B^G. It is even possible to use exactly the same description of the polarizability as in the full-potential calculations (see sections 2 and 4): using the diagonal term of equation (7) for $x\gt {{x}_{C}}={{R}_{C}}/\sigma$, replacing ${{{\rm cos} }^{2}}\theta$ by its approximate ℓ-independent mean value 1/2, and keeping a constant polarizability for $x\leqslant {{x}_{C}}$, we obtain (in reduced units)

$$\begin{eqnarray}&&{{C}^{G}}=-x_{00}^{4}/12+3x_{C}^{4}/48.\end{eqnarray} \tag{ 13 }$$

In order to determine bound levels and resonances, the calculation is performed in two steps. First, equation (9) is solved numerically by inward integration starting from a large value ${{x}_{\infty }}$, imposing only large x boundary conditions. For $\mathcal{E}\lt 0$, i.e., levels below threshold, this value has simply to be larger than the outer Condon point. For $\mathcal{E}\geqslant 0$, ${{x}_{\infty }}$ is chosen in the x-domain where the diagonal elements of the matrix ${\bf M}+\mathcal{E}{\mathbb{1}}$ reach their asymptotic form, that is $[\mathcal{E}-\ell (\ell +1)/{{x}^{2}}]$ for $\mathcal{E}\gt 0$ and $[1/{{x}^{6}}-\ell (\ell +1)/{{x}^{2}}]$ for $\mathcal{E}=0$. One can then use analytical solutions as initial values for the inward integration of the radial functions $y_{\ell }^{j}(x)$ in each channel and construct a set of linearly independent solutions ${{{\bf y}}^{j}}$ with the correct asymptotic behavior. There are n such solutions for bound levels and Siegert states. Their asymptotic behavior corresponds to either an exponentially decreasing function or an outgoing complex wave function in a given channel and zero in all others. For $\mathcal{E}\geqslant 0$, there are $2n$ linearly independent solutions, with an asymptotic behavior given by either a regular or an irregular Bessel function in a given channel and zero in all others, see reference [33] for $\epsilon$>0 and reference [35] for $\epsilon$=0.

In a second step, we calculate the physical solutions ${{{\bf z}}^{k}}$. They are linear combinations of the particular solutions ${{{\bf y}}^{j}}$ previously calculated. The coefficients are determined by imposing the radial components in each channel to vanish at the corresponding node position ${{x}_{0\ell }}$. This short range condition leads to a quantization of energy for the bound levels and Siegert states, see appendices A.1 and A.3, respectively. It also allows to determine the scattering length in the presence of the non-resonant field [36].

The continuous, n-fold degenerate spectrum at an energy $\mathcal{E}$ is described by multichannel scattering theory [37–39]. The chosen asymptotic boundary conditions allow for a direct determination of the energy-dependent reaction matrix, ${\bf K}(\mathcal{E})$, from which the scattering matrix, ${\bf S}(\mathcal{E})$, and the time-delay matrix, ${\bf Q}(\mathcal{E})$, are easily deduced. The details are presented in appendix A. In particular, the ${\bf Q}$ matrix is well-adapted to analyze shape resonances, by studying the energy variation of its lowest eigenvalue ${{q}_{1}}(\mathcal{E})$ which corresponds to Lorentzian profiles, see appendix A.2 for details. The eigenvalues ${\rm tan} [{{\tau }_{j}}]$ of the ${\bf K}$ matrix allow for determining the eigenphase sum $\tau (\mathcal{E})$. The energy variation of the derivative of the eigenphase sum yields also a profile of the shape resonances. The resonances are also finally characterized by calculating the energy variation of either population or mean value of $1/{{x}^{2}}$ inside the rotational barrier, see appendix A.2 for details.

To test the validity of the asymptotic model, we solve the asymptotic Schrödinger equation (7) and compare to results obtained with the full Hamiltonian (1) [17], using the ground state potential energy curve from [23], adjusted to yield the relevant scattering length. The polarizabilities are computed from equation (2) for $R\gt {{R}_{C}}=10$ a₀ with an atomic polarizability of ${{\alpha }_{0}}=186.25$a₀³ [42]; for $R\leqslant {{R}_{C}}$ the polarizabilities are taken to be constant. We first need to determine the nodal lines. To this end we use [13] which gives the energies and widths of shape resonances as a function of the position of a node at short range. Reversely, knowing the position of a field-free shape resonance, it is possible to find a node position ${{x}_{0\ell }}$ (in a chosen x-interval) that yields a resonance at this energy value. Starting from the positions of the field-free shape resonances ℓ = 4, 8, 12 and 16 [17], we first test nodal lines of the type (11). Since the coefficient A in equation (11) plays a minor role, it is taken to be constant and equal to $A={{A}^{G}}=-{{({{x}_{00}})}^{7}}/8$, the value of the 'universal' model [13]. $B(\ell )$ is taken to be a polynomial of degree 3 in $\ell (\ell +1)$; x₀₀ and $B(\ell )$ are determined by a fit to the field-free shape resonance positions (the degree of the polynomial is 3 to fit the 4 data points exactly). Note that this fit provides the correct value, $a=-2{{a}_{0}}$, of the field-free scattering length.

However, when using the ansatz (11) to determine, in addition to the field-free positions of the shape resonances, the slopes of their dependence on the non-resonant field, the result is disappointing: for the four resonances the slopes are smaller by a factor of approximately 1.75 compared to those obtained from the full Hamiltonian. This finding suggests that the contribution of the short-range part of the interaction with the non-resonant field (for $x\lt {{x}_{0\ell }}$) is non-negligible, rendering the use of field-independent nodal lines insufficient. Remarkably, the effect of the coupling at short range on the intensity dependence of the resonance positions can be simply compensated, at least roughly, by introducing a scaling factor in the field intensity.

The influence of the interaction with the non-resonant field at short range on the resonance positions can be fully accounted for in the asymptotic model by making the nodal lines intensity-dependent, cf equation (12). Assuming A to be constant, $A={{A}^{G}}$, as above, x₀₀ and $B(\ell )$, taken to be a polynomial of degree 4 in $\ell (\ell +1)$, are adjusted to reproduce exactly the nodes of the field-free wave functions with $\ell =0,4,8,12,16$. As above, this fit provides the correct value, $a=-2{{a}_{0}}$, of the field-free scattering length. Additionally, C, taken to be a polynomial of degree 3 in $\ell (\ell +1)$, is adjusted to reproduce exactly the variation of the node positions for $\ell =4,\;8,\;12,\;16$ with intensity when $\mathcal{I}$ is increased from 0 to 1 reduced unit. To this end, the ℓ-wave function of the rovibrational level closest to threshold is obtained numerically for $\mathcal{I}=0$ and $\mathcal{I}=1$, employing the Fourier grid method to solve equation (1) in a single channel approximation. For $\mathcal{I}=1$, the single channel calculation represents an approximation. It is, however, well justified by the very small ℓ mixing observed in a coupled channels calculation for $\mathcal{I}=1$, corresponding to $I=0.64\;$ GW cm⁻² for ⁸⁸Sr₂. The variation of the node positions at threshold with $\ell (\ell +1)$ is shown in figure 1 for $\mathcal{I}=0$. Also plotted are the node positions at threshold of the 'universal' model, i.e., $B(\ell )={{B}^{G}}=x_{00}^{5}/4$, corresponding to the same value of x₀₀.They do not deviate much from the node positions obtained from the full potential, except for large ℓ values. The node positions in the presence of a weak non-resonant field, $\mathcal{I}=1$ in reduced units, differ from those for $\mathcal{I}=0$ by only about $-0.5\times {{10}^{-4}}$ reduced units, not visible on the scale of the figure. Adding this small and simple linear intensity dependence to the nodal lines yields spectacular agreement of the asymptotic model with the full Hamiltonian. This is demonstrated by the upper panel of figure 2 which compares the results of the asymptotic model with intensity-dependent nodal lines to those of the full Hamiltonian: almost no difference is visible on the scale of the figure. A linear intensity dependence of the nodal lines thus allows for utilizing the asymptotic model up to very large field intensities.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Note that all crossings between resonances or levels in figure 2 are in fact avoided crossings, and the diabatized lines are simply labeled by the $\tilde{\ell }$ value equal to the field-free ℓ value. Figure 2 also shows the behavior of the resonance width as a function of field intensity for the example of $\tilde{\ell }=8$ (lower panel). The calculations using the asymptotic nodal line technique were performed with 11 coupled channels, but we have checked for several values of $\tilde{\ell }$ and $\mathcal{I}$ that the positions of the shape resonances (up to $\tilde{\ell }=20$) do not change when ℓ_max is increased (up to ${{\ell }_{{\rm max} }}=24$, corresponding to 13 coupled channels). The resonance positions and widths for $\tilde{\ell }=8,12,16$ shown in figure 2 have been obtained with the complex energy method, cf appendix A.3. For the resonance $\tilde{\ell }$ = 4, which is close to the top of the corresponding barrier at $\mathcal{I}=0$, the complex energy method does not apply and resonance profiles have been determined from the smallest eigenvalue of the time-delay matrix $Q(\mathcal{E})$, appendix A.2.

The intensity dependence of the resonance positions and widths shown in figure 2 is related to a strong hybridization of the rovibrational motion [16, 17]. The hybridization involves different aspects, which can be analyzed from profile calculations, cf appendix A.2. This is shown in figure 3, illustrating ℓ-mixing for the example of the $\tilde{\ell }=8$ resonance. The population density (per energy unit) trapped behind the centrifugal barrier (lower panel of figure 3) is essentially always concentrated in the $\ell =8$ channel. It increases rapidly when approaching the threshold. The crossing with the $\tilde{\ell }=12$ resonance does not visibly affect this evolution, whereas the crossing with $\tilde{\ell }=4$ involves a clear decrease of the population density in the $\ell =8$ channel (dip in the green line near 11 GW cm⁻²). Note the different behavior of the percentages at short and long range. In the short-range region the population is essentially concentrated in the $\ell =8$ channel; and the short-range percentages converge exactly to the population percentages in the different channels of the corresponding $\tilde{\ell }=8$ bound level at threshold (middle panel of figure 3). In contrast, the asymptotic percentages (upper panel of figure 3) are very different from the short-range ones, with a very small contribution of the $\ell =8$ partial wave and large contributions of partial waves with $\ell =0$ and 2 at high intensities. The asymptotic percentages represent the partial wave decomposition of the continuum wave function associated to the lowest eigenvalue of the time delay matrix $Q(\mathcal{E})$. We stress here that this wave function is, inside the multiply degenerate continuum, the only wave function exhibiting resonant behavior. The behavior of the asymptotic percentages would probably be important if dynamical processes were considered.

Zoom In Zoom Out Reset image size

The crucial free parameter in the asymptotic model, and the only free parameter in the 'universal' asymptotic model, is the value of the s-wave scattering length (in reduced units) which determines the node positions. It is thus particularly instructive to compare the ⁸⁸Sr₂ and ⁸⁶Sr⁸⁸Sr dimers. Since for small ℓ the differences between the 'universal' and the 'realistic' nodal lines are small, see figure 1, we use here 'universal' nodal lines (12) for ⁸⁶Sr⁸⁸Sr, with coefficients $A={{A}^{G}}$, $B={{B}^{G}}$, C = 0 and the value of x₀₀ determined by the s-wave scattering length. In this essentially explorative work, we have also limited the number of channels to 5 (ℓ_max = 8), sufficient to study the $\tilde{\ell }=2$ and $\tilde{\ell }=6$ resonances.

Encouraged by the very good agreement between the asymptotic model and the full Hamiltonian for the shape resonances, we calculate for both isotopomers, in addition to the shape resonances, bound levels very close to threshold. Figure 4 displays the positions of shape resonances and bound levels for ⁸⁸Sr₂, whereas the corresponding results for ⁸⁶Sr⁸⁸Sr are shown in figure 5. In both figures, we characterize each level or resonance by a value $\tilde{\ell }$, with the labeling done by continuity through avoided crossings ('diabatic' labeling) in two concording ways. First, we observe how the levels appear, as the number of channels in the calculations is enlarged; second, we analyze the channel decomposition for $\mathcal{I}=0$. The intensity dependence of the resonance and bound level positions is extremely different for ⁸⁸Sr₂ and ⁸⁶Sr⁸⁸Sr: in figure 5, new $\tilde{\ell }$ values, $\tilde{\ell }=2$ and $\tilde{\ell }=6$, appear, all crossings are widely avoided and the $\tilde{\ell }=6$ resonance crosses twice the threshold. It is worth mentioning that the theoretical energies obtained for the field-free bound ⁸⁸Sr₂ levels with $\ell =0$ and $\ell =2$ (−74.64 reduced units or −134.4 MHz and −36.37 reduced units or −65.5 MHz) are in good agreement with the experimental values of −136.7 and −66.6 MHz for the v = 62, $\ell =0$ and the v = 62, $\ell =2$ shape resonances, respectively [24]. In the presence of the non-resonant field, the bound levels in figures 4 and 5 are of three types: (i)'pure' bound levels which are bound also in the field-free case; (ii) bound levels which appear when a shape resonance is pushed below threshold as the non-resonant field intensity is increased; and (iii)'supplementary' bound levels, which start tangentially to the threshold, i.e., regular scattering states that become bound as the field intensity is increased. The latter are due to a deepening of the $\tilde{\ell }=0$ adiabatic potential as $\mathcal{I}$ increases. This effect is also observed for a very strong static electric field coupling to a permanent dipole moment [43].

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

Figure 6 analyzes the hybridization of a shape resonance for the example of the $\tilde{\ell }=6$ resonance in ⁸⁶Sr⁸⁸Sr. The features are similar to those shown in figure 3: in particular a drastic increase of the population density, especially in the $\ell =6$ channel, is observed when the resonance comes close to the threshold. Simultaneously, there is almost no more contribution of the $\ell =6$ partial wave in the asymptotic behavior.

Zoom In Zoom Out Reset image size

For a given value of energy, $\mathcal{E}=-{{k}^{2}}$, n linearly-independent solutions ${{{\bf y}}^{j}}$ are obtained by inward integration. Each solution j is related to a specific channel ${{\ell }_{j}}=2(j-1)$ by imposing $y_{\ell }^{j}(x)$ to behave asymptotically as $y_{\ell }^{j}(x)\propto {\rm exp} (-kx)$ in the ℓ_j channel and zero in all others. The physical solution ${\bf z}$ is a linear combination,

$$\begin{eqnarray}&&{\bf z}=\mathop{\sum }\limits_{j=1}^{n}{{a}_{j}}{{{\bf y}}^{j}},\end{eqnarray} \tag{ A1 }$$

where the radial component in each channel ℓ must vanish on the corresponding nodal line ${{x}_{0\ell }}$. The resulting linear system of n equations with n unknown variables a_j has a non-trivial solution if and only if

$$\begin{eqnarray}&&{{D}_{{\rm bound}}}(\mathcal{E})={\rm det} \left( y_{\ell }^{j}({{x}_{0\ell }}) \right)=0,\end{eqnarray} \tag{ A2 }$$

where ${\rm det} ({\bf M})$ is the determinant of the matrix ${\bf M}$. Equation (A2) is solved either by iteration on the energy $\mathcal{E}$ or by interpolation of ${{D}_{{\rm bound}}}(\mathcal{E})$ on a set of $\mathcal{E}$-values and finding the corresponding zeros, ${{\mathcal{E}}_{{\rm bound}}}$. Solution of the linear system of the n equations $\sum _{j=1}^{n}{{a}_{j}}\;y_{\ell }^{j}({{x}_{0\ell }})=0$ corresponding to the n ℓ-values at an energy ${{\mathcal{E}}_{{\rm bound}}}$ yields the coefficients a_j and thus the bound state wave functions of the various ℓ-channels. The coupled wave function ${\bf z}$ at the energy ${{\mathcal{E}}_{{\rm bound}}}$ is normalized to one, such that hybridization can be measured by the weights ${{\varpi }_{\ell }}({{\mathcal{E}}_{{\rm bound}}})=\int _{{{x}_{0\ell }}}^{{{x}_{\infty }}}{{[{{\sum }_{j}}{{a}_{j}}y_{\ell }^{j}(x)]}^{2}}\;{\rm d}x$.

To analyze the profiles of shape resonances in the n-fold degenerate continuous spectrum, we use multichannel scattering theory [37, 38]. For each energy $\mathcal{E}$, $\mathcal{E}={{k}^{2}}\gt 0$, we calculate $2n$ particular, linearly-independent, energy-normalized solutions of the Schrödinger equation (9) by inward integration. The initial conditions are taken at large distance ${{x}_{\infty }}$, where the centrifugal term $1/{{x}^{2}}$ prevails. For each channel ℓ_j (where again ${{\ell }_{j}}=2(j-1)$), we determine two particular solutions, denoted by ${{{\bf j}}^{j}}(x)$ and ${{{\bf y}}^{j}}(x)$, respectively, by imposing as asymptotic behavior in this channel either a regular $\sqrt{(k\pi x)/2}{{J}_{\ell +1/2}}(kx)$ or an irregular $\sqrt{(k\pi x)/2}{{Y}_{\ell +1/2}}(kx)$ energy-normalized Bessel function and zero in all other channels. The physical solutions of equation (9) are n linear combinations of the $2n$ calculated particular solutions which vanish on the nodal lines in each channel. Among all the possible sets of n particular combinations, we choose the 'standard' ones, ${{{\bf z}}^{j}}$, which asymptotically contain a regular component in the channel ℓ_j only

$$\begin{eqnarray}&&{{{\bf z}}^{j}}(x)=\mathop{\sum }\limits_{j^{\prime} =1}^{n}[\;{{\delta }_{j^{\prime} ,j}}\;{{{\bf j}}^{{{j}^{\prime }}}}(x)+\;{\bf K}_{j^{\prime} }^{j}\;{{{\bf y}}^{{{j}^{\prime }}}}(x)],\end{eqnarray} \tag{ A3 }$$

where ${\bf K}$ is the so-called reaction matrix [37, 38]. Introducing two n × n matrices

$$\begin{eqnarray}&&({{{\bf M}}_{{\rm reg}}})_{\ell }^{j}=j_{\ell }^{j}({{x}_{0\ell }}),\end{eqnarray} \tag{ A4a }$$

$$\begin{eqnarray}&&({{{\bf M}}_{{\rm irreg}}})_{\ell }^{j}=y_{\ell }^{j}({{x}_{0\ell }}),\end{eqnarray} \tag{ A4b }$$

the conditions that the wave functions vanish on the nodal lines allow us to determine the ${\bf K}-$ and ${\bf S}-$ matrices [37, 38]

$$\begin{eqnarray}&&{\bf K}=-{{({{{\bf M}}_{{\rm irreg}}})}^{-1}}\cdot {{{\bf M}}_{{\rm reg}}},\end{eqnarray} \tag{ A5 }$$

$$\begin{eqnarray}&&{\bf S}=({\bf 1}+\imath {\bf K})\cdot {{({\bf 1}-\imath {\bf K})}^{-1}}.\end{eqnarray} \tag{ A6 }$$

The existence and properties of a shape resonance can be determined by several different methods. An example of four different profiles that we have obtained in two particular cases (⁸⁶Sr⁸⁸Sr, $\tilde{\ell }=2$, at a field intensity $\mathcal{I}=5$ reduced units and ⁸⁸Sr₂, $\tilde{\ell }=8$, at a field intensity $\mathcal{I}=6.5$ reduced units) is displayed in figure A1 .

Zoom In Zoom Out Reset image size

The most suitable method for characterizing resonances is surely to study the time-delay matrix ${\bf Q}(\mathcal{E})$ [38, 39], which is related to the scattering matrix ${\bf S}$ by

$$\begin{eqnarray}&&{\bf Q}=-\imath \;{{{\bf S}}^{\dagger }}\;\cdot \;\frac{{\rm d}S}{{\rm d}\mathcal{E}}.\end{eqnarray} \tag{ A7 }$$

When a narrow and isolated resonance shape resonance is present, the energy dependence of the lowest eigenvalue of ${\bf Q}$, ${{q}_{1}}(\mathcal{E})$ exhibits a Lorentzian profile

$$\begin{eqnarray}&&{{q}_{1}}(\mathcal{E})\sim -\frac{{{\gamma }_{r}}}{{{(\mathcal{E}-{{\mathcal{E}}_{r}})}^{2}}+{{(\frac{1}{2}{{\gamma }_{r}})}^{2}}},\end{eqnarray} \tag{ A8 }$$

where ${{\mathcal{E}}_{r}}$ is the resonance energy and ${{\gamma }_{r}}$ its FWHM. The wavepacket associated to the corresponding eigenvector is resonantly delayed during its scattering by the attractive potential. This method can be extended to the analysis of overlapping resonances [47].

The lifetime $T=\hbar /\epsilon t$ of the resonance (in SI units) is calculated from the reduced lifetime $t=1/{{\gamma }_{r}}$ and the reduced unit of time $\hbar /\epsilon$. The channel-mixing of the resonance can be characterized by the eigenvector corresponding to ${{q}_{1}}({{\mathcal{E}}_{r}})$, which gives the partial wave decomposition of the continuum wave function (one among the n wave functions of the $n-$ multiple continuum) which concentrates the resonant character of the scattering at this intensity value.

The second method consists in diagonalizing the ${\bf K}$-matrix, with eigenvalues ${\rm tan} (\pi {{\tau }_{j}})$. The corresponding eigenvectors can be used to construct the so-called 'eigenchannel wave functions', ${{{\bf u}}^{j}}$, wich have the same asymptotic behavior, ${{u}^{j}}(x)\propto {{J}_{\ell +1/2}}(kx)+{\rm tan} (\pi {{\tau }_{j}}){{Y}_{\ell +1/2}}(kx)$, in all channels. The ${{\tau }_{j}}$ are called the eigenphase shifts. The total eigenphase shift

$$\begin{eqnarray}&&\tau (\mathcal{E})=\mathop{\sum }\limits_{j=1}^{n}\;{{\tau }_{j}},\end{eqnarray} \tag{ A9 }$$

increases by π when $\mathcal{E}$ passes through the resonance energy. The derivative with respect to the energy of the total eigenphase shift, ${{\tau }^{\prime }}(\mathcal{E})$, exhibits a resonance profile, since it is related to the trace of the ${\bf Q}$ matrix by ${\rm Tr}[{\bf Q}(\mathcal{E})]=2{{\tau }^{\prime }}(\mathcal{E})$.

It is finally possible to characterize the profile of a shape resonance from the radial components $v_{\ell }^{j}(x)$ of any orthonormalized set of continuum wave functions. In practice, we have used energy-normalized wave functions associated to eigenvalues of the matrix ${{{\bf K}}^{T}}\cdot {\bf K}$. Introduction of the transposed matrix ${{{\bf K}}^{T}}$ allows one to eliminate numerical problems related to small asymmetries of the matrix ${\bf K}$. In addition, the scalar product of the standard functions ${{{\bf z}}^{j}}$ defined in equation (A3) is equal to ${\mathbb{1}}+{{{\bf K}}^{T}}\cdot {\bf K}$. Since a shape resonance is a metastable state in which two atoms are temporally kept close to each other, a resonance profile is also expected for the density (per energy unit) inside the barrier and for the expectation value of the $1/{{x}^{2}}$. Precisely we calculate the following integrals

$$\begin{eqnarray}&&{{I}_{p}}=\mathop{\sum }\limits_{j=1}^{n}\;\mathop{\sum }\limits_{\ell =1}^{n}\int _{{{x}_{0\ell }}}^{x_{{\rm max} }^{\ell }}\;\frac{{{[z_{\ell }^{j}(x)]}^{2}}}{{{x}^{p}}}{\rm d}x.\end{eqnarray} \tag{ A10 }$$

With either p = 0 or p = 2; $x_{{\rm max} }^{\ell }$ is a very large value for p = 2; for p = 0, it is the position of the top of the centrifugal barrier $x_{{\rm top}}^{\ell }={{[\ell (\ell +1)/3]}^{-1/4}}$ for $\ell \gt 0$ and it is taken as $x_{{\rm top}}^{\ell =2}$ for $\ell =0$.

As shown in figure A1, the four calculated profiles, ${{q}_{1}}(\mathcal{E})$, ${{\tau }^{\prime }}(\mathcal{E})$, ${{I}_{2}}(\mathcal{E})$ and ${{I}_{0}}(\mathcal{E})$, exhibit similar shape, especially for narrow resonances. In the latter case, the profiles ${{q}_{1}}(\mathcal{E})$ are perfectly given by equation (A8).

Shape resonances with large ℓ are very narrow, even for low field intensity. This is due to the presence of the broad and high potential barrier. Resonances lying very close to the field-shifted dissociation limit also have a very small width. It is quite difficult to detect narrow resonances and to calculate their characteristics from an analysis of the resonance profiles as described in appendix A.2. As an alternative, we therefore calculate the resonances as Siegert states with complex energy ${{\mathcal{E}}_{S}}=\mathcal{E}-\imath \gamma /2$ where the real and imaginary part are related to the resonance energy and width [45, 46]. Siegert states are described by a complex wave function ${\bf z}$ whose asymptotic behavior corresponds to an outgoing wave in each channel.

To determine the Siegert states, we proceed similarly as for bound levels, cf appendix A.1. We first determine n linearly-independent particular complex solutions ${{{\bf y}}^{j}}$ of equation (9). Changing $\mathcal{E}$ into ${{\mathcal{E}}_{S}}={{k}^{2}}$ results in a complex k-value. Inward integration, imposing an outgoing wave asymptotically in the ℓ_j channel and zero in all others yields the ${{{\bf y}}^{j}}$. The outgoing wave is written as a combination of the regular and irregular Bessel functions, ${{(-1)}^{(\ell +1)/2}}\sqrt{\pi x/2}\left[ {{J}_{\ell +1/2}}(kx)+\imath {{Y}_{\ell +1/2}}(kx) \right]$. The physical Siegert wave function, obtained as a linear combination of the ${{{\bf y}}^{j}}$ solutions, has to satisfy the boundary condition at small x, i.e., the radial components in all channels have to vanish on the corresponding nodal line ${{x}_{0\ell }}$, calculated at the energy $\operatorname{Re}({{\mathcal{E}}_{S}})$. These conditions are equivalent to a vanishing determinant ${{D}_{{\rm siegert}}}({{\mathcal{E}}_{S}})$ of the radial ℓ-components of the n particular solutions $y_{\ell }^{j}(x)$ at the node positions. This condition quantifies the resonance energy to the value ${{\mathcal{E}}_{S,r}}={{\mathcal{E}}_{r}}-\imath {{\gamma }_{r}}/2$.