Dual Kinetic Balance Approach to the Dirac Equation for Axially Symmetric Systems: Application to Static and Time-Dependent Fields
E.B. Rozenbaum, D.A. Glazov, V.M. Shabaev, K.E. Sosnova, D.A. Telnov
DDual Kinetic Balance Approach to the Dirac Equation for Axially SymmetricSystems: Application to Static and Time-Dependent Fields
E.B. Rozenbaum, ∗ D.A. Glazov, V.M. Shabaev, K.E. Sosnova, and D.A. Telnov
Department of Physics, St. Petersburg State University,
St. Petersburg, Russia (Dated: October 9, 2018)Dual kinetic balance (DKB) technique was previously developed to eliminate spurious states inthe finite-basis-set-based solution of the Dirac equation in central fields. In the present paper, itis extended to the Dirac equation for systems with axial symmetry. The efficiency of the methodis demonstrated by the calculation of the energy spectra of hydrogenlike ions in presence of staticuniform electric or magnetic fields. In addition, the DKB basis set is implemented to solve the time-dependent Dirac equation making use of the split-operator technique. The excitation and ionizationprobabilities for the hydrogenlike argon and tin ions exposed to laser pulses are evaluated.
PACS numbers: 31.15.-p, 31.30.J-, 32.80.-t, 32.60.+i
INTRODUCTION
The finite-basis-set methods are widely used in atomic,molecular and solid state physics. These methods gen-erally possess high level of the numerical efficiency thatincludes a fast growth of the accuracy with increasingnumber of basis functions. It is known, however, that thestraightforward application of this kind of methods to theDirac equation leads to appearance of so-called spuriousstates (see Refs. [1–10] and references therein). As it wasshown in Refs. [7, 8], the spurious states originate fromthe restriction of the basis set to a finite number of func-tions. Several methods have been developed to solve thisproblem in the case of spherical symmetry. Among themis the well-known kinetic balance method [11–15] that im-plies the construction of the lower-component basis func-tions by applying the non-relativistic limit of the radialDirac operator to the upper component basis functions.In Refs. [16, 17], the basis was composed of the Gaus-sian spinors satisfying the boundary conditions for thecase of the finite nucleus. The Gaussian spinors obey thekinetic-balance condition (except for the non-relativisticlimit). Advantage of these approaches is a high accu-racy in calculating the bound-state energies. However,since the kinetic balance method violates the symmetryin treatment of the positive- and negative-energy states,application of this technique can be rather problematic ifthe contribution of the negative-energy continuum is sig-nificant. In particular, it takes place in calculations of theQED effects (e.g. for precise calculations of the g -factor).An effective and easily implementable method to get ridof the spurious states keeping the symmetry between theelectron and the positron states was proposed in Ref. [7].This method is called the dual kinetic balance (DKB) ap-proach. The efficiency of the DKB method was provedby the calculations of various relativistic and QED effectsin atomic systems [18–27]. In the present work, the DKBmethod is generalized to the Dirac equation for systemswith axial symmetry.At first, the finite-basis-set method is developed for stationary Dirac equation with axially symmetric poten-tial. The basis set is constructed from the one-componentbasis functions of radial and angular variables and trans-formed to the DKB form. To test the procedure, thespectra of hydrogenlike ions are obtained with makingno use of spherical symmetry of these ions. It is shownthat the DKB approach allows one to get rid of the spu-rious states while retaining the proper energy spectrum.To demonstrate the efficiency of the method, the spec-tra of hydrogenlike ions in strong uniform electric andmagnetic fields are calculated. The DKB approach doeseliminate the spurious solutions in these cases as well.The results for the Zeeman and Stark shifts of the levelswith n = 1 for Z = 1 , , , and 92 are compared withthose of the independent relativistic calculations basedon the perturbation theory.The developed method has a wide range of possibleapplications. In particular, it can be used for calcula-tions of the Zeeman effect, including linear ( g factor) andnon-linear contributions in magnetic fields. The latterappears to be important for middle- Z boronlike ions atthe present level of experimental accuracy [28]. In thisrespect, the DKB method represents a competitive alter-native to the traditionally employed perturbation theory.The time-dependent problems have been drawing muchmore attention during last years due to the rapid develop-ment of the laser technologies. There are several state-of-the-art laser facilities operating nowadays (see, e.g., [29])that provide extremely high intensities or frequencies ofthe radiation. Thus the processes involving the strong-field ionization and excitation are of a great interest [30–35]. Highly charged ions are among the most interestingobjects that can be experimentally studied with theselasers. Theoretical treatment of highly charged ions ex-posed to strong laser fields requires the fully relativisticconsideration. Thus the time-dependent Dirac equationis to be solved. Within the mostly relevant dipole ap-proximation and for the linearly polarized laser fields,these problems possess the axial symmetry and the so-lution can be based on the approach developed in the a r X i v : . [ phy s i c s . a t o m - ph ] N ov present paper for the stationary Dirac equation. Theparticular scheme of the solution of the time-dependentDirac equation is based on the split-operator techniqueand requires some transformations of the matrices. Thedeveloped technique is applied to calculate the excita-tion probabilities for hydrogenlike argon ion and the ion-ization probabilities for hydrogenlike tin ion exposed tovarious intense laser pulses. Some of the obtained re-sults for the excitation probabilities are compared withthose of the independent calculation (based on the first-order non-stationary perturbation theory). The resultsfor the ionization probability dependence on the laserwavelength are compared with the corresponding datafrom Ref. [35].Throughout the paper we assume (cid:126) = 1. STATIONARY DIRAC EQUATION
We consider the stationary Dirac equation H Ψ( r ) = E Ψ( r ) , (1)where H = c ( α · p ) + mc β + V , (2) V = V nucl ( r ) + V ext ( r, θ ) . (3)Here V nucl ( r ) is the nuclear potential. The stationaryexternal field potential V ext ( r, θ ), being dependent on theradial and angular variables ( r and θ ), takes differentforms for electric and magnetic fields. In the case ofexternal electric field, it is given in the length gauge by: V ext = − ( E · d ) = − e ( E · r ) , (4)where e is the electron charge and d is the dipole momentoperator: d = e r . For external magnetic field, we have V ext = − e ( α · A ) , (5)which can be rewritten in the particular gauge, A =12 [ H × r ], as: V ext = − e [r × α ] · H ) . (6)Here and below the external field is assumed to be di-rected along the z axis: E , H (cid:20) e z ; r , θ , and ϕ arethe corresponding spherical coordinates. In the case ofan axially symmetric field V ( r, θ ), the total angular mo-mentum J is not conserved. At the same time, the z -projection of the total angular momentum m J is con-served because the corresponding operator J z commuteswith the Hamiltonian:[ J z , H ] = 0 . (7) Consequently, H and J z have a common set of eigen-functions with explicit dependence on the azimuthal an-gle ϕ and thus the Dirac four-component wave function(bispinor) can be represented in the spherical coordinatesas follows: Ψ( r ) = 1 r G ( r, θ )e i ( m J − ) ϕ G ( r, θ )e i ( m J + ) ϕ ˙ ıF ( r, θ )e i ( m J − ) ϕ ˙ ıF ( r, θ )e i ( m J + ) ϕ . (8)Substitution of the form (8) into the Dirac equation (1)yields the equation H m J Φ = E Φ (9)for the function Φ( r, θ ) = G ( r, θ ) G ( r, θ ) F ( r, θ ) F ( r, θ ) . (10)In the case of external electric field, the operator H m J takes the following form: H m J = (cid:18) mc + V cD m J − cD m J − mc + V (cid:19) , (11)where V is given by Eqs. (3) and (4), D m J = ( σ z cos θ + σ x sin θ ) (cid:18) ∂∂r − r (cid:19) + 1 r ( σ x cos θ − σ z sin θ ) ∂∂θ + 1 r sin θ (cid:18) im J σ y + 12 σ x (cid:19) , (12) σ x , σ y , and σ z are the Pauli matrices.In the case of external magnetic field, the Hamiltonian H m J has the following form: H m J = mc + V nucl c (cid:16) D m J + ˜ D (cid:17) − c (cid:16) D m J + ˜ D (cid:17) − mc + V nucl , (13)where ˜ D = − e c H r sin θ ˙ ıσ y . (14)We note that D m J and ˜ D are anti-Hermitian opera-tors: D † m J = − D m J , (15)˜ D † = − ˜ D . (16)The scalar product in the space of the functions Φ isdefined by (cid:104) Φ a | Φ b (cid:105) = ∞ (cid:90) dr π (cid:90) dθ sin θ ( G a G b + G a G b + F a F b + F a F b ) . (17)Setting the boundary conditions,Φ( r, θ ) | r =0 = lim r →∞ Φ( r, θ ) = 0 , (18)leads Eq. (9) to be equivalent to the variational principle δ S = 0 for the functional S = (cid:68) Φ | H m J | Φ (cid:69) − E (cid:104) Φ | Φ (cid:105) , (19)where the undefined Lagrange factor E has the physicalmeaning of energy.Implementation of any kind of methods based on fi-nite basis sets starts with an approximate representationof unknown function Φ as a finite linear combination ofthe basis functions. Let N be the number of the four-component basis functions depending on the radial andangular variables ( r and θ ). We introduce a set of func-tions { W i ( r, θ ) }| Ni =1 , where r ∈ [0 , r max ] and θ ∈ [0 , π ].Then the function Φ can be expanded as follows:Φ( r, θ ) ∼ = N (cid:88) i =1 C i W i ( r, θ ) , (20)where C i are the expansion coefficients.By substitution of the expansion (20) into the varia-tional principle δ S = 0, the latter can be represented asa set of algebraic equations for the coefficients C i : d S dC i = 0 . (21)This system leads to the following generalized eigenvalueproblem: H ij C j = ES ij C j , (22)where the summation over the repeated indices is im-plied, H ij = ∞ (cid:90) dr π (cid:90) dθ sin θ [ W i ( r, θ )] † H m J W j ( r, θ ) , (23) S ij = ∞ (cid:90) dr π (cid:90) dθ sin θ [ W i ( r, θ )] † W j ( r, θ ) , (24) and the Hamiltonian H m J is defined by Eq. (11) or byEq. (13).Consider the construction of the basis set. Let N r and N θ be the numbers of the one-component basisfunctions depending on the r and θ variables, respec-tively. We denote these sets of functions as { B i r ( r ) } N r i r =1 and { Q i θ ( θ ) } N θ i θ =1 . The indices i r = 1 , . . . , N r , i θ =1 , . . . , N θ , and u = 1 , . . . , i = 1 , . . . , N introduced before ( N = 4 N r N θ ) as follows: i = ( u − N r N θ + ( i θ − N θ + i r . (25)Using these one-component single-variable function sets,we can construct the set of four-component functions W i ( r, θ ) = W ( u ) i r i θ ( r, θ ) of two variables. Then the ex-pansion (20) will take the form:Φ( r, θ ) ∼ = (cid:88) u =1 N r (cid:88) i r =1 N θ (cid:88) i θ =1 C ui r i θ W ( u ) i r i θ ( r, θ ) , (26)and indices i and j in Eqs. (21) - (24) should be replacedwith { i r , i θ , u } and { j r , j θ , v } , respectively.A straightforward way to construct the four-component basis functions depending on two variables( r and θ ) is: W ( u ) i r i θ ( r, θ ) = B i r ( r ) Q i θ ( θ ) e u , (27)wheree = , e = , e = , e = . (28)Our calculations with various finite-basis-set techniques,including the B-splines-based spectral approach [4–6] andthe generalized pseudo-spectral method [36], show thatthe basis (27) leads to the appearance of the spuriousstates.Following the idea of the DKB method, we should im-pose specific relations between the upper and lower com-ponents of the Dirac bispinor. These relations are derivedfrom the non-relativistic limit of the Dirac equation andgive the following basis functions in the case of axial sym-metry: W ( u ) i r i θ ( r, θ ) = Λ B i r ( r ) Q i θ ( θ ) e u , u = 1 , . . . , , (29)where Λ = − mc D m J − mc D m J . (30)It should be noted that, as in the case of central fields[7], the DKB approach for axially symmetric systems canbe used for the extended charge nucleus only. The point-like nucleus case can be accessed by the extrapolation ofthe extended-nucleus results to vanishing nuclear size.The discussion above is given for arbitrary basis sets { B i r ( r ) } N r i r =1 and { Q i θ ( θ ) } N θ i θ =1 . In the present work,the particular choice of the one-component basis func-tions is made as follows. The B -splines of some order k form the set of the one-component r -dependent ba-sis functions, { B i r ( r ) } N r i r =1 . The Legendre polynomials (cid:26) P l (cid:18) π θ − (cid:19)(cid:27) N θ − l =0 of degrees l = 0 , . . . N θ − θ -dependent basis functions, sothat in the previous notations Q i θ ( θ ) ≡ P i θ − (cid:18) π θ − (cid:19) .To demonstrate the absence of the spurious states incalculations based on the DKB method, in Table I wepresent the energy spectrum of hydrogenlike tin ion ( Z =50), evaluated for the extended nucleus case (the modelof the uniformly charged sphere is employed) with theplain basis set (27) and with the DKB basis set (29). Thecalculations are performed for the projection of the totalangular momentum m J = − / n -thorder are computed from the energies and the wave func-tions of the ( n − αZ in any order of the field strength. Tables Vand VI show that the present DKB approach fully re-produces the perturbation theory results. In this case,however, one should keep in mind that, strictly speak-ing, there are no discrete energy levels for atom in uni-form electric field. Instead, we have the quasi-stationarystates. It happens due to the tunneling effect for initiallylocalized electron state.The basis set of 78 radial B-splines of order k = 9and 17 Legendre polynomials (of orders from 0 to 16) isenough to obtain all the results presented in Tables I − VI.
TABLE I. Energy spectrum (in r.u.) of H-like tin ion( Z = 50, R nucl = 4 .
655 fm).n DKB off DKB on Exact values1 0.93106324090 0.93106324090 0.931063240860.970722241162 0.98261372423 0.98261372423 0.982613724230.98261424969 0.98261424969 0.982614249690.98321813638 0.98321813638 0.983218136380.986596701133 0.99234087351 0.99234087351 0.992340873510.99234102938 0.99234102938 0.992341029370.99252042806 0.99252042806 0.992520428060.99252042806 0.99252042806 0.992520428060.99257642386 0.99257642386 0.992576423860.99302522647TABLE II. Binding energy (in a.u.) of the ground state( m J = − /
2) of hydrogen atom in uniform magnetic field H = 0 . ≈ . · T). For comparison, the valueobtained in Refs. [37, 38] (they both coincide to all the pre-sented digits) is given. The complete perturbation-theory re-sult ∆ E PT and the individual contributions are listed as well.The terms missing in the breakdown (the odd orders > >
12) are zero to all the presented digits.This work Refs. [37, 38] PT order PT − . − . − . − . . . − . . − . . − . TABLE III. Binding energies (in a.u.) of the ground ( m J = − /
2) state and the lowest m J = − / H . H , a.u. m J This work Refs. [37–39]1 − / − . − . abc − / − . − . ac − / − . − . ab − / − . − . ab Taken from: a Ref. [37]; b Ref. [38]; c Ref. [39].TABLE IV. Zeeman shifts (in a.u.) of the ground states ofH-like ions. ∆ E PT are the perturbation-theory values. Thecontributions of the orders > Z = 18 R nucl = 3 .
427 fm H = 6 · T m J ∆ E DKB
PT order ∆ E PT +1 / . . − / − . − . ± / ± . ± / . ± / ∓ . ± / − . Z = 50 R nucl = 4 .
655 fm H = 6 · T m J ∆ E DKB
PT order ∆ E PT +1 / . . − / − . − . ± / ± . ± / . ± / ∓ . ± / − . Z = 92 R nucl = 5 . H = 6 · T m J ∆ E DKB
PT order ∆ E PT +1 / . . − / − . − . ± / ± . ± / . ± / ∓ . ± / − . TIME-DEPENDENT DIRAC EQUATION
We consider the time-dependent Dirac equation:˙ ı ∂∂t Ψ( r , t ) = H ( t )Ψ( r , t ) , (31)where H ( t ) = H + V ( t ) , (32) H = c ( α · p ) + mc β + V nucl ( r ) , (33)and V ( r , t ) describes the interaction with an externaltime-dependent field. In the following, we restrict our TABLE V. Stark shifts (in a.u.) of the ground state of hy-drogen atom and H-like argon ion. ∆ E PT and ∆ E nrPT arethe relativistic and non-relativistic perturbation-theory val-ues. The contributions beyond the shown ones are zero to allthe presented digits. Z = 1 E , V/m ∆ E DKB · PT order ∆ E PT · ∆ E nrPT · · − . − . − . · − . − . − . − . − . − . − . · − . − . − . − . − . − . − . Z = 18 R nucl = 3 .
427 fm E , V/m ∆ E DKB · PT order ∆ E PT · ∆ E nrPT · · − . − . − . · − . − . − . · − . − . − . n = 2 energy levelsof H-like argon ion ( Z = 18, R nucl = 3 .
427 fm). ∆ E PTD are calculated according to the approximate formulas fromRef. [42] derived within perturbation theory for degeneratelevels. E , V/m Level ∆ E DKB ∆ E PTD s . . · p / − . − . p / . . s . . · p / − . − . p / . . consideration to the time-dependent electric field withinthe dipole approximation: V ( r , t ) = − ( F ( t ) · d ) , (34)where F is the strength of the external electric field and d is the operator of the dipole moment: d = e r . Weassume F to be linearly polarized along the z axis: F ( t ) = F ( t ) e z . (35)Let ∆ t be a small time step. Given the initial wavefunction Ψ( r , r , t + ∆ t ) ≈ exp ( − ˙ ı ∆ tH [ r , t + ∆ t/ r , t ) . (36)For the function Φ defined in the previous section,this equation can be written asΦ( r, θ ; t + ∆ t ) ≈ exp (cid:16) − ˙ ı ∆ tH m J [ r, θ ; t + ∆ t/ (cid:17) Φ( r, θ ; t ) . (37)The direct application of these equations within thefinite-basis-set approach would be extremely time con-suming. For this reason, one needs to use special meth-ods to reduce the efforts. We use the split-operator tech-nique [43]. The implementation of this technique in theframework of the finite-basis-set method described aboverequires, however, some modifications that are presentedbelow.The split-operator method consists in the propagatorfactorization, e.g., as follows:Ψ( r , t + ∆ t ) ≈ exp (cid:20) − ˙ ı ∆ t H (cid:21) exp (cid:18) − ˙ ı ∆ t V (cid:20) r , t + ∆ t (cid:21)(cid:19) × exp (cid:20) − ˙ ı ∆ t H (cid:21) Ψ( r , t ) . (38)The exponential of the unperturbed Hamiltonian H istime-independent, and thus can be calculated only onceby the spectral expansion:exp (cid:20) − ˙ ı ∆ t H (cid:21) = (cid:88) k exp (cid:18) − ˙ ı ∆ t E k (cid:19) | Ψ k (cid:105) (cid:104) Ψ k | , (39)where H Ψ k = E k Ψ k . (40)In order to calculate the spectral expansion (39), we in-troduce the matrix and eigenvectors: H L = S − / H S − / , (cid:126)C L = S / (cid:126)C , (41)so that, instead of the generalized eigenvalue problem(22), we get the ordinary one: H L (cid:126)C L = E (cid:126)C L . (42)In order to get the highest possible efficiency, the time-dependent part V (cid:2) r , t + ∆ t (cid:3) should be represented by adiagonal matrix. According to Eqs. (34) and (35), wecan represent the matrix V as V ( t ) = F ( t ) · V , (43)where the matrix elements of V are given by V uvi r i θ j r j θ = ∞ (cid:90) dr π (cid:90) dθ sin θ (cid:104) W ( u ) i r i θ ( r, θ ) (cid:105) † × ( r cos θ ) W ( v ) j r j θ ( r, θ ) . (44)Let us consider the eigenvalue problem for the matrix V L = S − / V S − / : V L (cid:126)υ k = u Lk (cid:126)υ k (45) and construct the matrix of the eigenvectors: υ = ( (cid:126)υ (cid:126)υ (cid:126)υ . . . (cid:126)υ N ) . (46)Since the matrix V L is Hermitian, the matrix υ is unitary( υ † = υ − ) and the matrix V LV = υ † V L υ (47)is diagonal. Let us also denote: H LV = υ † H L υ, (cid:126)C LV = υ † (cid:126)C L . (48)With H and S substituted from equations (23) and(24), respectively, the time-dependent Dirac equation(31) takes the form:˙ ıS ddt (cid:126)C ( t ) = ( H + F ( t ) · V ) (cid:126)C ( t ) . (49)Multiplying Eq. (49) by υ † S − / , we get˙ ı ddt υ † S / (cid:126)C ( t ) = υ † S − / ( H + F ( t ) · V ) × S − / υυ † S / (cid:126)C ( t ) (50)or, using the notations (47) and (48),˙ ı ddt (cid:126)C LV ( t ) = (cid:0) H LV + F ( t ) · V LV (cid:1) (cid:126)C LV ( t ) . (51)This equation is suitable for the split-operator method.The short-term propagation can be performed as (cid:126)C LV ( t + ∆ t ) = exp (cid:20) − ˙ ı ∆ t H LV (cid:21) × exp (cid:20) − ˙ ı ∆ t F (cid:18) t + ∆ t (cid:19) V LV (cid:21) × exp (cid:20) − ˙ ı ∆ t H LV (cid:21) (cid:126)C LV ( t ) , (52)where the exponential of H LV is obtained by the spectralexpansion:exp (cid:20) − ˙ ı ∆ t H LV (cid:21) = (cid:88) k exp (cid:18) − ˙ ı ∆ t E k (cid:19) (cid:126)C LVk (cid:16) (cid:126)C
LVk (cid:17) † (53)and the matrix (cid:18) exp (cid:20) − ˙ ı ∆ t F (cid:18) t + ∆ t (cid:19) V LV (cid:21)(cid:19) ij = δ ij exp (cid:20) − ˙ ı ∆ t F (cid:18) t + ∆ t (cid:19) u Li (cid:21) (54)is diagonal (see Eqs. (45) - (47)).In order to calculate the transition and ionization prob-abilities, we have to project the propagated state onto the P u l s e s p e c t r u m , a . u . Energy, keV → → → → → → FIG. 1. (Color online) The energy spectrum of the Gaussian-shaped laser pulse used in our calculations. Vertical sticksindicate the photon energies necessary for one-photon (green)and two-photon (blue) transitions. vector or the subspace of interest. For instance, to cal-culate the survival probability in the initial state Ψ i , wehave to calculate the scalar product: |(cid:104) Ψ( t ) | Ψ i (cid:105)| = (cid:12)(cid:12)(cid:12) (cid:126)C † ( t ) · S (cid:126)C i (cid:12)(cid:12)(cid:12) = (cid:12)(cid:12)(cid:12)(cid:12) (cid:126)C † ( t ) · (cid:16) S / (cid:17) † υυ † S / (cid:126)C i (cid:12)(cid:12)(cid:12)(cid:12) = (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) υ † S / (cid:126)C ( t ) (cid:17) † · (cid:16) υ † S / (cid:126)C i (cid:17)(cid:12)(cid:12)(cid:12)(cid:12) = (cid:12)(cid:12)(cid:12)(cid:12)(cid:16) (cid:126)C LV ( t ) (cid:17) † · (cid:126)C LVi (cid:12)(cid:12)(cid:12)(cid:12) . (55)The developed methods have been applied to solve tworepresentative problems. First, the transition probabili-ties in the hydrogenlike argon ion ( Z = 18) exposed toa short Gaussian-shaped laser pulse are calculated. Theform of the pulse is given by the following: F ( t ) = e z F exp (cid:18) − t τ (cid:19) sin ( ωt ) , (56)where ω = 4 . − , τ = 0 .
63 as, and the peak inten-sity is I = 6 . · W/cm . Fig. 1 shows the energyspectrum (i.e. the Fourier transform) of this pulse. Wenote that the spectrum is broad enough, so that we getnot only one-photon transitions, but two-photon ones aswell. In Fig. 2, we present the transition probabilitiesfrom the ground 1 s state to the excited states due tothe interaction with the laser pulse. For comparison, thecorresponding results of the first-order time-dependentperturbation theory are shown in the same figure. Theinitial 1 s state survival probability is P s = 0 . Z = 50) exposed T r a n s i t i o n p r o b a b ili t y -8 -7 -6 -5 -4 -3 -8 -7 -6 -5 -4 -3 Excited states s p / p / s p / p / d / d / s p / p / d / d / f / f / Numerical calculationFirst-order time-dependent perturbation theory
FIG. 2. Electric-dipole transition probabilities from the 1 s state to the excited states for one-electron argon ion exposedto a Gaussian-shaped laser pulse. I o n i z a t i o n p r o b a b ili t y -8 -7 -6 -5 -4 -3 -3 -8 -7 -6 -5 -4 -3 -3 Wavelength, nm
FIG. 3. (Color online) Total ionization probability for a tinion as a function of the laser pulse wavelength. Unshadedred squares connected with the dashed red line: the resultsof our calculations; shaded black points connected with thesolid black line: data taken from Ref. [35] to a short sin -shaped laser pulse. In this calculation, thepulse is chosen in the form: F ( t ) = e z F sin (cid:18) πtT (cid:19) sin ( ωt ) , t ∈ [0 , T ] , (57)where the wavelength λ and the pulse duration T can beexpressed through the carrier frequency as λ = 2 πc/ω and T = 2 πN /ω . The calculations are performed for N = 20 and the peak intensity I = 5 × W/cm .Fig. 3 displays the full ionization probability as a func-tion of the laser wavelength for all the other parametersof the system kept constant. This plot is in a good agree-ment with the corresponding data from Ref. [35]. CONCLUSION
The efficient and easily implementable DKB approachsolves the problem of the spurious states related to theuse of the finite basis sets for the Dirac equation. In thepresent paper, this method is generalized for the case ofthe axial symmetry. Generalized DKB method proved tobe accurate and stable in this case. It opens the new wayfor the fully relativistic theoretical treatment of both sta-tionary and time-dependent axially symmetric problems,e.g., of ions and atoms exposed to external fields. Theefficiency of the method is demonstrated by calculatingthe energies of hydrogenlike ions with non-perturbativeaccount for static uniform external electric or magneticfields. The Zeeman and Stark energy shifts are comparedwith the perturbation theory calculations. It is shownthat the higher orders of the perturbation theory expan-sion can be reproduced by the methods developed in thepresent paper.For the purpose of solving the time-dependent prob-lem, the finite basis set technique (not regarding theparticular choice of the basis set) was adapted to takeadvantage of the split-operator method by the transfor-mation of the matrix of the external potential into thediagonal representation. With this technique, the transi-tion and ionization probabilities for the ions exposed tothe laser pulses are evaluated. The results are comparedwith the corresponding data from other papers or withthe independently obtained values. The solution of thetime-dependent Dirac equation with the set of the dis-cussed approaches is shown to be correct, accurate andnumerically efficient.
ACKNOWLEDGEMENT
We thank A. S´aenz for providing us with the numericaldata from Ref. [35]. Valuable discussions with A. S´aenzand J. L. Sanz-Vicario are gratefully acknowledged. Thiswork was supported in part by RFBR (Grants No. 12-02-31803 and 13-02-00630), by the Ministry of Educationand Science of Russian Federation (Project No. 8420), bythe FAIR-Russia Research Center, and by the non-profit“Dynasty” Foundation. ∗ [email protected][1] V. F. Bratzev, G. B. Deyneka, and I. I. Tupitsyn, Izv.Acad. Nauk SSSR (in Russian), Ser. Fiz. , 2655 (1977)[Bull. Acad. Sci. USSR: Phys. Ser. , 173 (1977)][2] G. W. F. Drake and S. P. Goldman, Phys. Rev. A ,2093 (1981)[3] W. Kutzelnigg, Int. J. Quantum Chem. , 107 (1984) [4] P. Schl¨uter, G. Soff, K -H. Wietschorke, and W. Greiner,Journal of Physics B: Atomic and Molecular Physics ,1685 (1985)[5] W. R. Johnson and J. Sapirstein, Phys. Rev. Lett. ,1126 (1986)[6] W. R. Johnson, S. A. Blundell, and J. Sapirstein, Phys.Rev. A , 307 (1988)[7] V. M. Shabaev, I. I. Tupitsyn, V. A. Yerokhin, G. Plu-nien, and G. Soff, Phys. Rev. Lett. , 130405 (2004)[8] I. I. Tupitsyn and V. M. Shabaev, Optics and Spec-troscopy , 183 (2008)[9] M. Lewin and E. Sere, Proceedings of the London math-ematical society , 864 (2010)[10] M. Lewin and E. Sere, arXiv:1306.5401 (2013)[11] Y. S. Lee and A. D. McLean, J. Chem. Phys. , 735(1982)[12] I. P. Grant, Phys. Rev. A , 1230 (1982)[13] K. G. Dyall, I. P. Grant, and S. Wilson, J. Phys. B ,L45 (1984)[14] R. E. Stanton and S. Havriliak, J. Chem. Phys. , 1910(1984)[15] H. M. Quiney, I. P. Grant, and S. Wilson, Phys. Scr. ,460 (1987)[16] Y. Ishikawa and H. M. Quiney, Int J. Quantum Chem.Quantum Biol. Symp. , 523 (1987)[17] Y. Ishikawa, K. Koc, and W. H. E. Schwarz, Chem. Phys. , 239 (1997)[18] V. M. Shabaev, K. Pachucki, I. I. Tupitsyn, and V. A.Yerokhin, Phys. Rev. Lett. , 213002 (2005)[19] A. V. Volotka, D. A. Glazov, V. M. Shabaev, I. I. Tupit-syn, and G. Plunien, Phys. Rev. Lett. , 033005 (2009)[20] K. Beloy and A. Derevianko, Computer Physics Commu-nications , 310 (2008)[21] D. A. Glazov, A. V. Volotka, V. M. Shabaev, I. I. Tupit-syn, and G. Plunien, Phys. Rev. A , 062112 (2010)[22] V. A. Yerokhin, Phys. Rev. A , 012507 (2011)[23] Q. M. Sun, W. J. Liu, and W. Kutzelnigg, TheoreticalChemistry Accounts , 423 (2011)[24] E. A. Chernovskaya, O. Yu. Andreev, and L. N. Lab-zowsky, Phys. Rev. A , 062515 (2011)[25] J. Zatorski, N. S. Oreshkina, C. H. Keitel, and Z. Har-man, Phys. Rev. Lett. , 063005 (2012)[26] S. R. McConnell, A. N. Artemyev, M. Mai, andA. Surzhykov, Phys. Rev. A , 052705 (2012)[27] A. Ong, J. C. Berengut, and V. V. Flambaum, arXiv:1310.2685 (2013)[28] D. von Lindenfels, M. Wiesel, D. A. Glazov, A. V.Volotka, M. M. Sokolov, V. M. Shabaev, G. Plunien,W. Quint, G. Birkl, A. Martin, and M. Vogel, Phys. Rev.A , 023412 (2013)[29] M. Vogel, W. Quint, G. G. Paulus, and Th. St¨ohlker,Nucl. Instr. Meth. B , 65 (2012)[30] A. A. Sorokin, S. V. Bobashev, T. Feigl, K. Tiedtke,H. Wabnitz, and M. Richter, Phys. Rev. Lett. , 213002(2007)[31] M. Richter, M. Ya. Amusia, S. V. Bobashev, T. Feigl,P. N. Juranic, M. Martins, A. A. Sorokin, and K. Tiedtke,Phys. Rev. Lett. , 163002 (2009)[32] S. Selsto, E. Lindroth, and J. Bengtsson, Phys. Rev. A , 043418 (2009)[33] M. S. Pindzola, J. A. Ludlow, and J. Colgan, Phys. Rev.A , 063431 (2010)[34] M. S. Pindzola, Sh. A. Abdel-Naby, F. Robicheaux, andJ. Colgan, Phys. Rev. A , 032701 (2012) [35] Y. V. Vanne and A. Saenz, Phys. Rev. A , 033411(2012)[36] D. Gottlieb and S. A. Orszag, Numerical analysis of spec-tral methods: theory and applications (Society for Indus-trial and Applied Mathematics, Philadelphia, 1977)[37] Z. Chen and S. P. Goldman, Phys. Rev. A , 1722 (1992)[38] A. Rutkowski and A. Poszwa, Physica Scripta , 484(2005)[39] H. Nakashima and H. Nakatsuji, The Astrophysical Jour-nal , 528 (2010) [40] L. B. Mendelson, Phys. Rev. , 90 (1968)[41] V. D. Ovsiannikov and V. G. Pal’chikov, Can. J. Phys. , 1401 (2002)[42] S. A. Zapryagaev, Optics and Spectroscopy , 892(1978)[43] J. A. Fleck, J. R. Morris, and M. D. Feit, Appl. phys.10