Recursive algorithm for arrays of generalized Bessel functions: Numerical access to Dirac-Volkov solutions

In the relativistic and the nonrelativistic theoretical treatment of moderate and high-power laser-matter interaction, the generalized Bessel function occurs naturally when a Schrödinger-Volkov and Dirac-Volkov solution is expanded into plane waves. For the evaluation of cross sections of quantum electrodynamic processes in a linearly polarized laser field, it is often necessary to evaluate large arrays of generalized Bessel functions, of arbitrary index but with fixed arguments. We show that the generalized Bessel function can be evaluated, in a numerically stable way, by utilizing a recurrence relation and a normalization condition only, without having to compute any initial value. We demonstrate the utility of the method by illustrating the quantum-classical correspondence of the Dirac-Volkov solutions via numerical calculations.

Figure 1

(Color online) Illustration of the different saddle-point regions of $J_n(x,y)$, for the two qualitatively different cases: number 1, with $8y>x$ (here $x=y={10}^3$ was used), and number 2, with $8y<x$ $(x=10y={10}^3)$. In case 1 [panel (a)], the transition from region $b_1$, where only one saddle point is real, to $c_1$, where two real saddle points contribute, occurs precisely at $n=-2y+x=-{10}^3$ (see also Table ). The complex oscillating behavior in region $c_1$ can be understood as interference between the contributions from the two real saddle points. In case 2 [panel (b)], we have $8y<x$, with only three as opposed to four qualitatively different regions (see also Table ).

Figure 2

(Color online) The five-term recurrence relation (11) has four linearly independent solutions. Note the logarithmic scale. The values $x=y={10}^3$ were used for the calculation, corresponding to case 1 $(8y>x)$. The solutions are labeled by $J_n$ [red line, the true generalized Bessel function $J_n(x,y)$], $Y_n$ (light blue line), $X_n$ (green line), and $Z_n$ (blue line). The numerically obtained solutions $X_n$, $Y_n$, and $Z_n$ have been shifted vertically by multiplication with an appropriate constant (of order ${10}^{10}$ for $Y_n$, ${10}^{-10}$ for $Z_n$, and ${10}^{-20}$ for $X_n$) for clarity. The separation of the different saddle-point regions is marked with dashed lines. The regions $a_1$, $b_1$, $c_1$, and $d_1$ are described in Table .

Figure 3

(Color online) Panel (a) shows the absolute value of the coefficients from Eqs. (27, 28), used for the transformed recurrence relations (25, 26). Here, a starting index of $M_{-}=-4000$ and initial values of ${\xi }_{M_{-}}^{1,2,3}={\lambda }_{M_{-}}^{1,2}=1$ were used. Note that for better visibility, all curves except ${\xi }_n^1$ have been shifted vertically on the logarithmic ordinate axis by multiplication with a suitable constant (${10}^{10}$ for ${\xi }_n^3$, ${10}^5$ for ${\xi }_n^2$, ${10}^{-7}$ for ${\lambda }_n^1$, and ${10}^{-12}$ for ${\lambda }_n^2$). In panel (b), we display the absolute values of $f_n$ and $g_n$, which result after completing step 3 in the algorithm described, i.e., before normalizing. The dash-dotted lines indicate the starting indices $M_{\pm }$, here $M_{-}=-3300$ and $M_{+}=2350$. The cutoff indices $n_{-}=-3000$, $n_{+}=2063$ are plotted with dashed lines. An example of a suitable “matching index” $K=0$ (see step 4 of the algorithm) is drawn by a solid line. The initial values used to calculate the curves were $f_{M_{+}}={10}^{-20}∕2$, $f_{M_{+}+1}=0$, $f_{M_{+}+2}={10}^{-20}$, and $g_{M_{+}}=0$, $g_{M_{+}+1}={10}^{17}$. Note that in panel (b), no vertical shifting was applied. The inset shows a magnification of the cutoff region for positive $n$, where the diverging behavior of $g_n$ for $n>n_{+}$ is clearly seen. In both graphs we have $x=y={10}^3$, same as in Fig. 2.

Figure 4

(Color online) The relative error ${ϵ}_{\mathrm{rel}}$, as defined in Eq. (43), as a function of the “safety margin” parameter $\Delta$ [see Eq. (42)]. In each of the different parameter ranges considered, an exponential decrease of the obtained error with the safety margin parameter is observed, demonstrating the applicability of the recursive method. In (a) we consider parameters such that $x=2y$, in (b) we have a large ratio $y/x$, whereas in (c), we have a small ratio $y/x$, and in (d), results for the cutoff region are presented. A detailed explanation of the parameter regions considered is in the text. In all parts, the black circles represent the approximation for the relative error obtained from Eqs. (45, 46).

Figure 5

(Color online) Illustration of the classical-quantum correspondence of a laser-dressed electron. Shown to the left [(a)] with a solid blue line is the classical energy $u^0$ [see Eq. (57)] as a function of the phase $\varphi$, for $\omega =1\mathrm{eV}$, laser intensity $I={10}^{16}\mathrm{W}∕{\mathrm{cm}}^2$ (corresponding to $\mid ea\mid ∕m=0.1$ where $a$ is the laser polarization four-vector), initial energy $p^0=2m$, and initial angle $\theta =0.54°$, with $\vec{p}\cdot \vec{k}=\mid \vec{p}\mid \omega \cos \theta$. These parameters give $x=3.3\times {10}^3$, and $y=-2.7\times {10}^3$ for the arguments of $J_n(x,y)$. The dashed lines show the minimum classical energy $u_{\min }^0$, the maximum classical energy $u_{\max }^0$, the average energy $q^0$, and the intermediate energy level (local maximum) $u_{\mathrm{int}}^0$, as indicated in the center of the figure. In (b), we display the quantum mechanical amplitude $J_n(x,y)$ of energy level $n$, which has energy $q^0+n\omega$ [see Eq. (49)]. Here positive index $n$ corresponds to absorbing $n$ number of photons from the laser field, while negative $n$ means emitting $\mid n\mid$ number of photons into the laser mode. The graph is arranged such that $n=0$ corresponds to the average energy $q^0$. The cutoff indices are nicely reproduced by the classical maxima and minima.