Contents lists available at ScienceDirect

Radiation Physics and Chemistry

journal homepage: www.elsevier.com/locate/radphyschem

S Radiation Physics and Chemistry

Theory and calculation of the atomic photoeffect

Lorenzo Sabbatuccia,b, Francesc Salvatb,n

a Department of Industrial Engineering (DIN), Laboratory of Montecuccolino, Alma Mater Studiorum University of Bologna, via dei Colli 16, 40136 Bologna, Italy

b Facultat de Física (ECM and ICC), Universitat de Barcelona, Diagonal 645, 08028 Barcelona, Spain

CrossMark

HIGHLIGHTS

• The theory of the atomic photoeffect is presented in a form suited for calculation.

• Explicit formulas are given for both ionization and excitation to bound levels.

• Simple calculation method for photoelectron angular distributions.

• Finite level widths can be included by convolution with a Lorentz profile.

ARTICLE INFO ABSTRACT

The so-called elementary theory of the atomic photoeffect is presented in a form that is suited for practical numerical calculation of subshell cross sections and angular distributions of emitted photo-electrons. Atomic states are described within the independent-electron approximation, with bound and free one-electron orbitals that are solutions of the Dirac equation with the Dirac-Hartree-Fock-Slater self-consistent potential of the ground-state configuration. Detailed derivations are given of subshell cross sections for both excitation to discrete bound levels and ionization. In the case of ionization, the cross section differential in the direction of the photoelectron is obtained for partially polarized photons, with the polarization specified by means of the Stokes parameters. The theoretical formulas have been implemented in a computer program named PHOTACS that calculates tables of excitation and ionization cross sections for any element and subshell. Numerical calculations are practicable for excitations to final states with the principal quantum number up to about 20 and for ionization by photons with energy up to about 2 MeV. Elaborate extrapolation schemes for determining the subshell cross section for excitation to bound levels with larger principal quantum numbers and for ionization by photons with higher energies are described. The effect of the finite width of atomic energy levels is accounted for by convolving the calculated subshell cross-section with a Lorentzian profile.

© 2015 The Authors. Published by Elsevier Ltd. This is an open access article under the CC BY-NC-ND

license (http://creativecommons.org/licenses/by-nc-nd/4.0/).

Available online 22 October 2015_

Keywords: Photoeffect

Subshell cross sections Angular distribution of photoelectrons Photon polarization Independent-electron approximation Dirac-Hartree-Fock-Slater self-consistent potential

1. Introduction

In the photoelectric effect, or photoeffect, a photon is absorbed by a target atom and, as a result, an atomic electron is emitted or promoted to a bound open orbital thus leaving the residual ion or atom in an excited state. The latter subsequently decays to its ground state through a cascade of radiative and non-radiative transitions with emission of characteristic X-rays and Auger electrons. For photons with intermediate and low energies, the pho-toeffect dominates the transfer of energy from the photon field to charged particles. The so-called elementary theory of the process

* Corresponding author. E-mail address: francesc.salvat@ub.edu (F. Salvat).

has been described by Pratt et al. (1973) and by Scofield (1973). In their formulation the states of the atom are approximated by a model of independent electrons in a common central potential, and the interaction between the target atom and the electromagnetic field is treated as a perturbation to first order. This approach neglects electron correlations, i.e., the collective character of the response of atomic electrons to the external field. Correlation effects have been studied by many authors using a variety of theoretical methods such as the random-phase approximation (Johnson and Lin, 1979), many-body perturbation theory (Chang and Fano, 1976), R-matrix perturbation theory (Burke and Taylor, 1975), and the time-dependent local-density approximation (Zangwill and Soven, 1980; Parpia and Johnson, 1983; Liberman and Zangwill, 1984). Calculations with these methods are quite

http://dx.doi.org/10.1016/j.radphyschem.2015.10.021

0969-806X/© 2015 The Authors. Published by Elsevier Ltd. This is an open access article under the CC BY-NC-ND license (http://creativecommons.org/licenses/by-nc-nd/4.0/).

involved, and numerical results have been published only for specific atoms and limited energy ranges.

Quantitative information on the photoeffect is required for practical applications (e.g., X-ray fluorescence, X-ray photoelectron spectroscopy), as well as for Monte Carlo simulation of photon transport. Numerical tables of subshell cross sections (for ioniza-tion and for excitation to bound levels) and atomic cross sections are included in the Evaluated Photon Data Library (EPDL) (Cullen et al., 1997). The XCOM program (Berger et al., 2005) gives atomic cross sections for photoionization essentially equivalent to those in the EPDL. Both the EPDL and XCOM databases are based on calculations performed by Scofield (1973) using the independent-electron model with the self-consistent Dirac-Hartree-Fock-Slater (DHFS) potential. They are considered to be the most reliable source of general information available to date; indeed, practically all modern Monte Carlo codes for photon transport utilize the EPDL. A systematic comparison of Scofield's cross sections with experimental data has been made by Saloman et al. (1988). It is worth mentioning that the theory, as well as the numerical tables calculated from it, applies to free atoms. Differences are to be expected for molecules and solids, particularly near absorption edges, partly because of aggregation effects on the atomic potential and also because of the EXAFS effect (extended X-ray absorption fine structure) (Lee et al., 1981).

Calculations based on approximate independent-electron models, such as the DHFS model, are affected by possible inaccuracies of the adopted central potential. A simple strategy to account for inaccuracies in the atomic potential is provided by the normalization screening approximation of Pratt and co-workers (Pratt, 1960a; Schmickley and Pratt, 1967; Pratt and Tseng, 1972; Pratt et al., 1973). According to this approximation, the subshell cross sections calculated from the DHFS potential and from a more elaborate atomic model (e.g.., the multi-configuration Dirac-Fock self-consistent model implemented in the program of Desclaux, 1975, 1977) differ essentially by a constant factor, which is equal to the ratio of electron densities near the nucleus. That is, in principle, one can improve the DHFS cross sections by multiplying them by an energy-independent factor, which is readily obtained from atomic-structure calculations that are more elaborate than the DHFS ones.

The angular distribution of photoelectrons is needed in X-ray photoelectron spectroscopy, (Hemmers et al., 2004; Fujikawa et al., 2007; Kovér, 2010), as well as in radiation transport calculations (Fernández et al., 1993; Salvat, 2015). Surprisingly, information on the angular distribution of photoelectrons, consistent with the subshell cross sections, is quite limited, or unavailable (see, e.g., Trzhaskovskaya et al., 2006 and references therein). As a matter of fact, most Monte Carlo photon transport codes still rely on the Sauter (1931) formula, which gives the differential cross section for the ground state of hydrogenic ions obtained from the plane-wave Born approximation. The most elaborate tables of angular distributions available for all subshells of the elements are those given by Trzhaskovskaya et al. (2001, 2002, 2006), which were calculated within the quadrupole approximation and parameterized using the formulas proposed by Cooper (1990, 1993). These tables apply to unpolarized and linearly polarized photons and cover the energy range from 100 eV to 5 keV. Theoretical studies of angular distributions for soft-X-ray absorption by Derevianko et al. (2000) and Amusia et al. (2001 ) revealed the importance of higher-mul-tipole corrections.

The present paper describes the elementary theory of the photoeffect in a concise but complete form, organized so as to allow systematic evaluations of cross sections for both ionization and excitation to bound levels. Following Scofield and others, we consider a model of independent electrons in the DHFS self-consistent potential. Although Scofield's numerical results are known

to be fairly accurate, the size of his database is moderate and, in some energy ranges, the grid energies are too widely spaced to ensure that interpolation errors are less than the numerical accuracy of the tabulated data. Nowadays, instead of interpolating from a limited database, subshell ionization cross sections can be calculated interactively, even on modest personal computers, using numerical algorithms that are highly accurate. We have written a robust, flexible Fortran program, named photacs, which calculates subshell cross sections for arbitrary atomic potentials. The calculations converge for excitation to levels with principal quantum number up to about 20 and for ionization by photons with energies up to about 2 MeV. Results from the program do confirm the accuracy of Scofield's data, but also reveal near-threshold features that are invisible in Scofield's database because of the large spacing of its energy grid. photacs includes tables of theoretical and empirical atomic level widths and allows the calculated subshell cross section (excitation plus ionization) to be transformed into a continuous function of the photon energy, through a convolution with a Lorentzian line profile. Once the photoemission cross section for a given photon energy is calculated, the angular distribution of photoelectrons can be obtained with little additional effort because radial integrals, which take most of the numerical work, were already computed to obtain the cross section. The code photacs optionally calculates the angular distribution of photoelectrons emitted as a result of the absorption of partially polarized photons in a given subshell. Arbitrary photon polarization is described by means of the density matrix expressed in terms of the stokes parameters.

The paper is organized as follows. section 2 presents an overview of the theory of the photoeffect. We describe a systematic extrapolation scheme to account for excitations to bound levels near the ionization threshold, as well as an analytical formula for extrapolating the photoionization cross section to arbitrarily high energies. Finally, the effect of the finite width of atomic energy levels and the normalization screening correction are considered. To lighten the Theory section, notation conventions and relevant mathematical formulas are summarized in Appendix A. In Section 3 we derive general formulas for the angular distribution of pho-toelectrons released by photons with arbitrary polarizations. photon polarization is described by means of the density matrix and the associated Stokes parameters, following the conventions presented in Appendix B. The practical calculation of cross sections and photoelectron angular distributions is considered in Section 4, which contains a brief presentation of the program photacs and results from several illustrative calculations. Finally, in Section 5 we provide some concluding comments.

2. Theory

We consider here the elementary theory of the photoelectric effect as described by Pratt et al. (1973) and by Scofield (1973). The states of the target atom are represented as single Slater determinants, built with one-electron orbitals tynKm (see Section A.1 in Appendix A) that are solutions of the Dirac equation for a central potential V(r). This potential is set equal to the self-consistent DHFS potential for the ground-state configuration (see Section A.2), which gives one-electron energy eigenvalues close to the experimental binding energies. As mentioned above, the theory assumes that the interaction between the atom and the electromagnetic field is weak and, therefore, it can be treated as a perturbation to first order. consequently, the theory is able to describe the process in which a single photon of energy W = hm is absorbed. Multiple photon absorption which occurs, for instance, in highly focused laser beams will not be considered.

Conventional time-dependent perturbation theory (Baym, 1974; Akhiezer and Berestetskii, 1965) leads to results equivalent to the one-active-electron approximation. That is, in calculations of the photoeffect we may assume that the active electron jumps from a bound orbital = yn K m of energy ea to a final orbital with energy eb = ea + W. The other (non-active) atomic electrons are treated as mere spectators; their only effect is to screen the nuclear charge, which is accounted for by the atomic potential V(r). The same potential is used to describe both the initial and final states and, consequently, all the orbitals involved in the calculations are mutually orthogonal. It is worth mentioning here that independent-electron models are not expected to yield accurate results for excitation to bound levels, because the adopted excitation energies (determined from energy eigenvalues of Dirac's one-electron wave equation) may differ largely from their actual values. Nonetheless, these models provide a fairly accurate description of photoionization, as well as useful insight into the gross features of the excitation spectrum of free atoms.

Let us start by considering a beam of photons with energy W = hm and wave vector k (k = m/c), in a pure polarization state Z, incident on a target atom at the origin of coordinates. The absorption of a photon causes the transition of the active electron from the initial orbital % with energy ea < 0, and ionization energy Ea = |eaI, to an excited orbital with energy eb = ea + W, which may be either bound (when W < Ea, excitation) or free (if W > Ea, ionization). To simplify calculations, we assume that photons propagate along the direction of the z-axis. This implies that the polarization vector lies on the x-y plane and, hence, the state Z can be represented as (see Appendix B)

Z = cos(a/2) + sin(a/2) exp(iy9)e2,

where the unit vectors e1 and e2 represent, respectively, states of linear polarization in the directions of the x and y axes. The cases ft = 0, n and ft = ± n/2 correspond to linear and circular polarizations, respectively. We would like to point out that the spherical components of the polarization vector [Z = S,VZ, see Eq. (A.35)] are such that & = 0 (transverse field) and |f-1|2 + |f+1|2 = 1 (normalization).

The cross section for excitation of the active electron to a bound orbita1, ¥b = Vnbw is given by

0 6x0,1 (W z)

vnbkbmb;nakama \vv ''s)

(2n)2e2ch |

\'vlnbkbmb;na kama

(W,Z)5( -ea - W),

with the transition matrix element

'v'nbkbmb;na kam

a (W, Z)= Z Gexc(Kb, mb; Ka, ma),

c(Kb, mb; Ka, ma ) = (Wnbkbmb \a exp (ikr)jl^naKama)

The formulas (A.40) in Appendix A allow the calculation of the spherical components of the vector Gexc as sums of products of radial integrals and matrix elements of Racah tensors [see also (A.44)].

The delta distribution in Eq. (2) indicates that absorption is possible only at resonance, i.e., when the excitation energy eb - ea coincides with the energy of the photon. A finite absorption cross section, which varies continuously with the photon energy, is obtained when the energy width of atomic levels is taken into account (see Section 2.5).

In the case of ionization, the final orbital belongs to the continuum spectrum ( eb > 0). The differential cross section (DCS) for absorption of a photon, with emission of the active electron with spin mSb, linear momentum fokb, and direction of motion within

the solid angle element dkb about the direction kb, is (Pratt et al., 1973; Scofield, 1973)

d _ion,1 dab;na Kama

(W, Z ) =

(2n)2e2 kb(eb + meC2)

MbOSaKama (W, Z)) dkb (5)

where the transition matrix element

b;na kama

(W, Z) = Z\1^»!a exP (ik r)^naKama)

is to be evaluated on the energy shell, i.e., with eb = ea + W. The result (5) is obtained from Fermi's golden rule (see, e.g., Baym, 1974), with the final state of the electron represented by a distorted plane wave (DPW), Eq. (A.15), with incoming spherical distortion, % = ^k-)msi, and normalized in the form (A.17), so that the density of states per unit volume in kb - space is unity. Inserting the expansion (A.15) of the DPW, the matrix element becomes

Ayi'on

b;na Kama

(W,Z ) = Til

eb + 2meC2 n (eb + meC2)

2 i-lb exp (i<5

X Z Gion(Kb, mb; Ka, ma

mSb °Kbmb

G ion(Kb, mb; Ka, ma ) = (yebkbmb \SC exp ( ikr)| y/nakama)

are matrix elements of the operator a exp(ik r) in the basis of spherical waves, analogous to those of excitation, Eq. (4). The spherical components of Gion can be calculated by means of formulas (A.44) in Appendix A.

2.1. Cross sections of closed subshells

In practice, target atoms are randomly oriented and the measurable quantity is the cross section for photoabsorption by the electrons of a subshell (naKa)qa with qa electrons. Calculations are easier in the case of closed subshells, with qa = 2 |Ka| electrons, because of various sum rules obeyed by the angular integrals in transition matrix elements. For open subshells, with qa < 2|Ka |, the average over atomic orientations is equivalent to assuming that the individual orbitals in the subshell have a fractional occupation number equal to qa/2|Ka|; the resulting partial cross section is then equal to that of the closed subshell times qa/2\Ka|.

In the case of excitation, the measured cross section includes transitions to any state of the final energy level eb. Generally, the final level is empty (qb=0), although occasionally it may correspond to an open subshell (nbKb)qb filled with qb ( 0 < qb < 2|Kb |) electrons. Hence, the cross section for one-electron excitation from the subshell (naKa)qa to a bound level eb by absorption of a photon of polarization Z is given by

tfgf (W ) =

Qa 21Kbl - qb 2|Ka| 2 |Kb|

aexU (w Z)

nbkbmb;na kama\vv ' b )

a (W ) 5 (eb -ea - W ),

where we have introduced the "reduced" cross section, 0a (W), defined by

(W) = (2n)2e2ch 2K"l - qbgexc(Kb, Ka),

W 2 |kJ 2\Kb\ v '

9exc(ib, Ka) - XX \MnHHbmb;naKama ma mb

W, Z ))

x ( - 1)vZvG-xc(Kb, mb; Ka, ma)

is the sum of squared transition-matrix elements over orbitals of the initial and final energy levels, which is given by Eq. (A.47),

9exc( Kb, Ka)- x ux^eaka] + [ ""j^cka]

The quantities e,m^nare defined by Eqs. (A.43) with radial

integrals m<RnbKb.eaKa evaluated for the final bound orbital VnbKbmb.

Note that the reduced cross section (10) is independent of the polarization of the absorbed photon. Naturally, this result is a direct consequence of the spherical symmetry of closed subshells.

Let us now consider the observed DCS for ionization of a subshell (naKa)qa. In addition to the assumption that the target atom is randomly oriented, we suppose that the final spin state of the emitted photoelectron is not measured. The observed DCS is thus obtained by summing over degenerate initial states and over final spin states, mSb,

d<n(W, Z) (2n)2e2 kb( + meC2) qa

21 Kal

xXX iMOUma(W,z)|2.

ma mSb

Inserting the expansion (7), we have d<n(W, Z) (2n)2e2ch kb qa

W n6b 2|Ka|

x x i x i-lb exp(i^Kj)Z Gion(Kb, mb; Ka, ma)

x i x i-lb exp(i^Kb)Z-Gion(Kb, mb; Ka, ma

k-b,mb

i XmSbXmSb

b (kb),

where the sum in square brackets is equal to the 2 x 2 unit matrix. The total cross section for ionization of the subshell (naKa)qa is obtained by integrating the DCS over the directions of the emerging photoelectron:

. , N rdo'r(W, Z) a <n(W)- f a ' dkb.

Using the orthogonality of the spherical spinors, Qi^mfo (kb), we obtain

^ (W J-SW* £ X (Kb, „ ),

9ion (Kb, Ka ) - xx |Z Gion ( Kb mb; Kama ))

x (- 1)vZ-vCT(Kb, mb; Ka

is the sum of squared transition-matrix elements over the orbitals

of the active subshell and over the spherical waves with the final energy, which is given by Eq. (A.47),

§^ (кb, Ka) - x i [e^6bib;caia]2 + [J^a]

with the quantities e,m^i

' nbKb;eaKa defined by Eqs. (A.43). We see that, as in the case of excitation, the total cross section is independent of the polarization of the absorbed photon.

(12) 2.2. The dipole approximation

The calculation of cross sections for both excitation and ioni-zation by photons with low energies, such that the wavelength X = 2n/k is much larger than the average radial distance of the electrons in the active subshell, can be simplified by using the dipole approximation. In this approximation, the quantity kr is assumed to be much less than unity, and the exponential in the matrix elements (4) and (8) is replaced by unity:

Gdip(ib, mb; Ka, ma) ^ (ye„ kbmb \Wnakam^j■

Using the commutation relation of the Dirac Hamiltonian, Eq. (A.1), with r,

[ «, r]

i ch a,

we get the familiar "length" form of the matrix elements in the dipole approximation:

Gdip(ib, mb; Ka, ma) - -h W (Web.Kb.mb\a,

l, ka, ma

Recalling that the spherical components of the vector r are rv = rC1v(r), the angular integrals of these matrix elements reduce to matrix elements of the Racah tensors. We have

Cdip (ib, mb; Ka, ma) - ch W {QKb,mb\Ci,v\QKa,m^ ®e with the radial integrals

r [ PebKb (r) PeaKa ( r) + QebKb ( r) QeaKa ( r)]r dr.

Within the dipole approximation, the quantities (12) and (18) reduce to

Qd'P(Kb, Ka) - xx

- 3(ch)2

x ( - 1)vZv C-Vp(Kb, mb; Ka, ma) -±1

lbyJb C(1) layJa ) [ oeb,kb-,ea,Ka] ■

Although the dipole approximation yields realistic values of the total cross section for long-wavelength photons (see Section 4), the DCS (i.e., the angular distribution of photoelectrons) obtained from this approximation is qualitatively correct only for photon energies near the ionization threshold (see Fig. 5).

2.3. Near-edge excitation cross section

Because of memory and time constraints, our computer program PHOTACS calculates reduced excitation cross sections, ^ (W), for bound levels with principal quantum number nb up to ncut = 18. The numerical value of the quantity §dip(Kb, Ka), Eq. (24), qualifies the levels e„bKb. Cross sections of dipole-forbidden levels [with §dip(Kb, Ka) = 0] have values that are typically several orders of magnitude smaller than those of dipole-allowed levels [with §dip(Kb, Ka) ^ 0]. To reduce the size of the tables, and to

ma k lm

simplify further calculations, photacs delivers cross-section values only for excitations to dipole-allowed levels. Cross sections for transitions to dipole-forbidden levels are calculated and added to the cross section of the nearest dipole-allowed level. In other words, the reduced cross section ff^swa (W) given by the program is the cross section of the allowed level enbKb plus small contributions of excitations to neighbouring forbidden levels.

For a given value of Kb, there exists an infinite series of bound levels e

nb kb,

whose energies increase monotonically with the principal quantum number nb and have an accumulation point at zero energy. To account for the cross sections for excitation to levels with nb larger than 18, we rely on the fact that, regarding the bound levels with energies close to zero as a quasi-continuum, the cross section for excitation continuously joins the ionization cross section at the edge, W = Ea (see, e.g., Watanabe, 1965; Fano and Cooper, 1968). Explicitly, from the cross sections for excitation to

the levels e

nb kb'

we define the following step function

Sa ( Kb; W)

;a (W )

enh enb-1

if enb-i - ea < W <enh - ea

where e„b are the midpoints of the energy intervals between successive energy levels:

enbkb + enb+1,kb

That is, the excitation cross section is spread over the energy interval that is closest to the final level. Numerical results confirm that the step function Sa (xb; W) does smoothly join the ionization cross section at the edge (see Fig. 3 below). The smooth matching of the discrete and the continuum is due to the fact that, for a given angular momentum Kb and for not too large radii, the radial functions of bound orbitals with very large nb differ from those of free orbitals with very small kinetic energies only by a normalization factor. When free-state orbitals are normalized on the energy scale

y^eVm' WeKm dr = S (e' - e)Sk'k 8„

and bound orbitals are normalized to unity, Eq. (A.11), the radial functions of bound and free states with energies close to zero are such that

b (r )

e nb e nb-1

for sufficiently small radii r, large nb, and small eb. This quantity, which coincides with the "spread factor" in the definition (25), can be regarded as the number of discrete levels per unit energy interval, the analogue to the density of states of the continuum.

In calculations of excitation to bound levels each series of (optically-allowed) levels e„bKb, with Kb fixed, is considered separately. Let n0 (Kb) denote the principal quantum number of the first level in the series. We describe the first three levels in the series [with nb = n0(Kb), n0(Kb)+ 1, n0(Kb) + 2] as discrete resonances. The levels with nb > n0 (Kb) + 2 are treated as a quasi-continuum, which extends from a certain cut-off energy W^ to W = Ea, with the excitation cross section expressed as a cubic polynomial in Ea - W:

ôaexc(kb; W) = A^,0 + Akb,l(Ea - W) + Akb,2(Ea - W)2 + AKb, 3 (Ea - W)3.

The coefficients of this polynomial are determined from a least-squares fit to a table of values of the step function (25) at the excitation energies of the levels, {enbKb, Sa (Kb; e%Kb - ea) } with

integral of the analytical approximation (29) over W in the interval from W£J to the upper limit of the subinterval of the nb=18 level equals the sum of excitation cross sections of the levels enbKb with nb = n0 (Kb) + 3 to 18. We assume that this polynomial represents the quasi-continuum faithfully and, hence, that it does account for excitations to levels with nb greater than 18. These considerations can be verified numerically using our computer programs. A graphical analysis of the case of excitations from the K shell of argon atoms is given below in Fig. 3.

Note that pseudo-continua of series with different Kb extend over different intervals. The cross section for excitations to bound levels is then given by

a™(W) = X OblKa(W) <5(eb - ea - W)

X â?c (Kb; W)(W - WKb)(Ea - W),

where 0 (x) is the unit step function (= 1 if x > 0, = 0 otherwise). The first summation is over discrete resonances (excitations to the three lowest levels of each Kb series), and the second summation accounts for excitations to levels in the pseudo-continua.

2.4. Ionization cross sections at high energies

The convergence rate of the series (16) decreases when the energy of the photon increases. Therefore, the program photacs is able to compute the ionization cross section ol°n (W) for photons with energy W up to a certain maximum value Waut. Cross sections for photons with energies higher than WCut can be obtained by extrapolating the calculated numerical cross sections using the approximate formula given by Pratt (1960a) (see also Pratt et al., 1973) for K-shell electrons, which combines Pratt's high-energy limit with the general energy dependence determined by Gavrila (1959). The formula for the cross section per electron in the K shell of an atom of atomic number Z reads

(W ) = âo

(fly)3 exp[- 2(a/y9)arccos a]

(W/mec2)4 a2

x M (fl){ 1 + na [ N (fl)/M (fl)] + R (a)}, ^

where a = Za [a = e2/(ch) ^ 1/137 is the fine-structure constant], | = 1 - V1 - a2, ¡3 is the velocity of the photoelectron in units of c

(W - Ea)(W - Ea + 2meC2)

(W - Ea + meC2)2

1 - fl2

W - Ea mec2 '

The quantity a0 is the high-energy cross section per electron, given

âo = 2na0 a5a3—-

(29) The functions M (fl) and N (fl) are defined by

m (fl) = 4 + Y^ v' 3 Y + 1

2yVY2 - 1 Y - -Jy2 - 1

Y + V Y2 - 1

N (fi)

63 25 8 - 4y + 34--+ + ^

Y Y2 Y3

15fi 3 15(y2 - 3Y + 2)

1 - fi

1 + fi

Finally, the quantity R(a) is defined so as to reproduce Pratt's high-energy limit of the total cross section (per electron) for the pho-toeffect in the K shell, which is given by Pratt (1960a) and Pratt et al. (1973).

lim ok (W)- oo F (a) with

F (a)-^ a_2i J /> /> ( ^ )"

_ -a-2«- , 2« 4 J-1

x [ a + i-,/1 - a2 x] x j (1 - x2y2)1-«

1-4+2«

- 2(1 - y2)(1 - x2y2)-

i ax +( a )2x2

Consistency of Eq. (31) with this high-energy limit implies that exp( - 2a cos-1 a)

F ( a)--

1--a ■

R ( a)

Values of the function F(a) calculated numerically by Pratt, with an estimated accuracy of 0.1 %, are given in Table I of Pratt (1960a). Numerical values of R(a) are also given by Pratt et al. (1973) in their Table 6.1. Unfortunately, the latter table contains several erroneous signs. To determine an accurate analytical approximation to R(a), we have calculated the function F(a) for a dense grid of a values by computing the double integral (38) using an adaptive Gauss-Le-gendre quadrature method, which allows control of numerical errors, to an accuracy of about 10-7. The calculated table extends from a=0.01-0.8, corresponding to Z up to 109, and agrees with Pratt's values (Pratt, 1960a) within the claimed accuracy of the latter. A least-squares fit of our tabulated values leads to the approximation

R(a) - 0.00372a - 0.16326a2 + 0.94375a3 - 0.71732a4,

which, when inserted into expression (39), approximates the numerical values of F(a) with relative error less than 0.01%.

As indicated by Pratt et al. (1973) the ratios of cross sections of different subshells are nearly energy independent. Hence, although the formula (31) was derived for the K shell, it can be employed to extrapolate the DHFS cross sections ofother subshells for energies above their numerical cut-off Wacut. In our calculations we set

1 (W )-

O Pratt (W^t )

oPratt (W )

for W > Wau In practice, this extrapolation scheme works very well for inner subshells with binding energies larger than about 100 eV, for which the extrapolation matches the numerical cross sections in a wide energy interval, and not so well for outer shells with Ea < 100 eV because the slopes of the numerical and the extrapolation curves at W - W,Cut are slightly different. However, the relative contributions of the outer subshells are less than about 10-5 and, consequently, the error introduced by the extrapolation has no practical effects. A similar extrapolation method was

employed by Hubbell et al. (1980), and also by Cullen et al. (1997), who used a semiempirical formula that reproduces Pratt's high-energy values for the K shell.

2.5. Finite level widths and experimental ionization energies

Up to this point we have been assuming that excited states are stationary. This assumption leads to the delusive delta function in the excitation cross section, Eq. (9). In reality, excited states decay by radiative and non-radiative transitions and, consequently, have a finite mean life Tb and a natural level width rb = h/rb. Assuming that excited estates decay exponentially with time, it is concluded that absorption lines have a Lorentz profile centred at the resonance energy, with full width at half-maximum equal to rb (see, e.g., Sakurai, 1967). Consequently, results from measurements of the energy E of an excited level eb are expected to follow the Lorentz distribution

L ((; 6b - E) - -

n (6b - E)2 + (b/2)2

In addition, the cross section for excitation of an electron from the level ea (corresponding to the ground state) to a level eb by absorption of photons with energy W should be calculated as the integral over the continuous level profile:

c (W)- f One£C6;a (W )S (E -6a - W)L (rb; 6b - E)dE

- OZb.a (W ) L (fb; 6b -6a - W ).

This result differs from the expression (9) for stationary levels only in that the Lorentz distribution replaces the delta resonance.

As noted by Richtmyer et al. (1934), in the case of excitations that produce a single vacancy in an inner subshell (naKa), the first stage of the subsequent decay of the atom is the filling of that vacancy by electrons from nearest subshells, a process practically independent of the condition of the excited electron. Consequently, all excited levels with a vacancy in subshell (na Ka) have approximately the same level width, ra, the so-called "core-level width". Calculated values of core-level widths of free atoms are given in the EADL (Perkins et al., 1991). Campbell and Papp (2001) give a complete set of recommended widths for K to N7 levels of atoms obtained from consideration of available experimental data. Typically, the core-level widths increase with the binding energy of the subshell, and are of the order of 0.1 eV or less for weakly bound subshells, and reach values of the order of 100 eV for K shells of transuranic elements.

Assuming that the level width is a characteristic of the active subshell (naKa)qa where the vacancy is created, its influence on the photoeffect can be accounted for by convolving the calculated cross section with the Lorentzian distribution (42). That is to say [cf. Eq. (30)],

Va (W) = x ^a (W) L (fa; 6b - 6a - W)

+ x fl oaexc (b; E)L (ra; E - W)dE

kbWcub

oaon(E)L(ra; E - W) dE,

where again, the first summation is over discrete resonances with excitation energies less than the corresponding cut-off Wacubt. The integrals of the second term on the right-hand side, with Vfc (Kb; W) expressed by the polynomial approximation (29), can be evaluated analytically. The last term can be evaluated similarly when the ionization cross section &aon (W) is obtained from a pre-calculated table of values, at suitably spaced energies, by means of

cubic spline interpolation.

There is a further empirical correction to be considered. Although the DHFS eigenvalues, Ea, are quite close to the experimental subshell ionization energies Efxp (Salvat and Fernández-Varea, 2009), the differences induce appreciable shifts of the absorption edges (typically, of a few eV). We can correct the calculated cross sections for this discrepancy, by simply shifting the energy scale and setting

af (W) = aa(W + Ea - Eaexp).

The post-processing program photacs-pp (see below) allows this energy shift to be applied automatically. The default experimental ionization energies Eexp have been taken from the compilation of Carlson (1975), which covers all the subshells of the elements from Z = 1-106. Alternatively, the program allows ionization energies to be used from the more updated compilation by Williams (2011), which is mostly based on X-ray photoelectron measurements on samples prepared in ultra-high vacuum conditions.

7 o 0.100

« S 0.010 [i]

: N 1-3 :

" N \\ \ \ ■i», \ \ \ \ \ -

= \ \ s \ Xj M5~" N7

" \ ; \ * \ \ \ L3 ' ___

\ \ > % \ VL1 MÍ M3 ----

: __ _____ ---:

: -

Fig. 1. Relative reduction of the subshell cross sections introduced by the normalization screening correction, obtained from electron density ratios resulting from MCDF and DHFS self-consistent atomic calculations.

2.6. Normalization screening approximation

A primary ingredient of the calculations is the DHFS potential, which provides only a rough approximation to the atomic wave functions. More elaborate atomic structure calculations employ the multi-configuration Dirac-Fock (MCDF) method (Desclaux, 1975). Unfortunately, the theoretical scheme underlying this method does not permit easy calculation of photoelectric cross sections, because the MCDF equations involve a non-local potential that is different for each subshell.

A simple method for estimating how the use of more accurate wave functions would affect the cross section for the photoeffect is provided by pratt's normalization screening approximation, which is presented by pratt (1960b) and discussed by Schmickley and Pratt (1967), Pratt and Tseng (1972), and Pratt et al. (1973). The key idea under this approximation is that, for photon energies a few keV above the ionization threshold, the dominant contributions to the transition matrix elements come from radial distances of the order of the reduced electron Compton wavelength, %e = h/mec = 3.8616 x 10-13m. In other words, from distances that are substantially larger than the radius of the nucleus but much smaller than the average radial distance of electrons in bound orbitals. At radii important for photoabsorption, the only effect of screening is a change in the normalization of bound states. It then follows that the photoelectric cross sections for the screened atomic potential (scr) and for the Coulomb potential of the bare nucleus (Coul) differ by an energy-independent factor

aacr(W) = 32 aC0UÍ(W),

where e2 is the ratio of electron densities at r=0. This ratio is determined by the normalization of the bound-state wave functions:

32 = lim

[P2aka (r) + Qiaka (^

[pn2aka (r) + qlaka (r)]c

lim r[P"aKa(r)]si

[P"aKa

( r )]c

The accuracy of the normalization screening approximation can be practically assessed by running our program with the DHFS potential and with the (unscreened) Coulomb potential. It should be noted, however, that the proportionality (46) only holds when a point nucleus is considered in the DHFS calculation. Otherwise, the radial functions for the DHFS potential and for the Coulomb potential have different shapes, a fact that complicates the theoretical analysis. Indeed, the approximation is found to work very well for photons with energies larger than the DHFS ionization energies. We may expect that it will perform even better for

"similar" screened potentials. Schmickley and Pratt (1967) went a step further and suggested that the normalization screening approximation can be effectively employed to relate cross sections calculated with different atomic models. Scofield (1973), more explicitly, stated that "the approximation probably holds for an atom as a whole, if the single-particle model is given up." Consistently, he calculates cross sections for all the elements with the DHFS potential and lists the normalization ratios e2 for orbitals obtained from the DHFS code and from restricted relativistic Hartree-Fock calculations for Z = 1-54.

Following Scofield, we have evaluated E2 for all subshells of the ground-state configuration of the elements Z = 1-99 using wave functions obtained from calculations with our DHFS program and with the MCDF program of Desclaux (1975, 1977). The resulting density ratios for the K shell and the L, M and N subshells of the elements are represented in Fig. 1. Because density ratios are less than, and close to unity, the displayed quantity is the difference 1 - EMMCDF,a/E|2HFS,a. We see that the effect of the normalization screening correction is a reduction of the cross section of the order of a few percent or less, except for outer subshells with small binding energies, where it can rise up to about 30 percent. Furthermore, the correction is larger for the outer subshells and, consequently, its effect is most visible for photons with relatively low energies. We can thus apply the normalization screening correction to the photoionization cross sections obtained from DHFS calculations

afr (W) =

3 DHFS,a

"MCDF,a aDHFS

to get better estimates of cross sections, whose accuracy is expected to be comparable to that of the MCDF model.

3. Angular distribution of photoelectrons

Although the polarization of the photon does not affect the total cross section, it does have an effect on the angular distribution of photoelectrons which, classically, accelerate in the direction of the electric field of the incident electromagnetic wave. The DCS for ionization, Eq. (14), can be expressed as

dâ™ (W, Z) _ (2n)2e2ch kb qa

W neb 2|Kal

X Z gKama (Z, kb) + gKama ( kb)

where we have introduced the emission amplitude (spinor)

SKama (^ k b) = z ^ eXP (iSKb )

X Z Gion(Kb, m6; Ka, ma))bmb (kb)

= z ukb,mb;ka,ma &kbmb (kb),

, = i-lb exp(iS-b)| Z ( - 1)v?vG-On(Kb, mb; -a, ma)|,

and [see Eq. (A.44a)]

Gil"(Kbmb; Kama) = ( - 1)ja-ma £ i;(j'a, jb, ma, - mt\j, + l)

X { ^JjfrKb;eaka ± i ^JbKb;eaka} (52)

with [see Eqs. (A.43b) and (A.43c)]

^Ebk b;eaka

ebkb'<eaka

iV4 j Mb je

ebk b;eaka'

The value of the component G0on is irrelevant here because the polarization vector is on the x-y plane, i.e., Z0 = 0.

Evidently, the DCS (49) is more difficult to compute than the total cross section, Eq. (16), because of the dependence of the emission amplitude on the magnetic quantum numbers of the initial and final states. As we will see, the DCS can be effectively computed from Eq. (49). Of course, the resulting DCS is consistent with the total cross (16), that is, when the emission amplitude is calculated using the same number of Kb terms as in expression (16), numerical evaluation of the integral (15) should yield the same result as Eq. (16).

To obtain a more explicit expansion of the DCS in terms of Le-gendre functions (similar, e.g., to that derived by Scofield, 1989 for linearly polarized photons), we introduce the expression (A.6) of the spherical spinors, and make use of the orthogonality of the unit spi-nors xH, to write

dâr(W, Z) _ (2n)2e2ch kb qa

neb 2|Ka|

ZZ(- 1)v

xX Z ilb-lbexp[i(SKb - SKb)]

11 KbKb

X Z G-v (k b, mb; K a, ma ) G* ' (Kb, mb; k a, ma ) mbmb

x (lb, 2-, mb - 1,1 jb, m^(f'b, 1 mb - 1, j'b, mbj mb-1 (k b)Y *b,mb-1 (k b|

fb,mb-H

Introducing the Clebsch-Gordan series for the spherical harmonics, and after rather lengthy transformations using conventional angular momentum recoupling algebra (Rose, 1995; Edmonds, 1960; Var-shalovich et al., 1988), the DCS can be written in the form (cf. Scofield, 1989)

dâ™(W , Z) _ (2n)2e2ch kb qa

neb 2 K

^z ZvZ*Z Al'Yi,v-v(kb)

a vv ' i ^ '

AV = (- 1)l+v', z (- 1 ) ja +jb +jb + 1/2

' 4n kbikbJ'

x ^(2jb + 1)(2j'b + 1)(^Vb + 1)(2lb + 1)(2i + 1)(2i' + 1)

'lb lb l^i jb jb ii i A

0 0 0J | lb lb 2II jb jb ja I

ilb-lb+- exp[i(SKb - SKb)] x { e y/ - i v m yJ }

x { e W + i v ' m yJ }

1 eb-b;£ak-a eb-b;eaka j '

v - v v - v

where (:::) and {:::} are, respectively, Wigner's 3j and 6j symbols (Rose, 1995; Edmonds, 1960; Varshalovich et al., 1988). A convenient simplification is obtained by using the formula

u (l, lb, lb)

l jb jb

lb lb l

V(2lb + 1)(2lb + 1)

'jb jb l ]

0 0 0J] lb lb 1

where u (l,4, lb) is the parity factor (A.30). Noting that ( - 1)v = - 1, and considering the symmetry of the 3j symbols, we get

Alv' = (- 1 ) ja +1/2

2l + 1 4n

x z M + 1)(2jb + 1)(2/ + 1)(2/ + 1)

kbjkjfi'

'l jb jb

0 - I 1

x u (l, lb, lb

x (J J' l |ilb-lb^W

I v - v ' v ' - v

i i' l|

jb jb ja I

exp[i(S-b - S-b

1 u eb-b;

- i v m VJ

eaka ^ eb-b;eaka

x { e W + i v ' m }

1 eb- b;£ak-a ^ebkfceakaj-

From the symmetry properties of the 3j and 6j symbols, it follows that

Al' = [ A^^v ]*. (58)

In addition, because of the parity factor in the reduced matrix element of the Racah tensor, Eq. (A.29), the vector coupling coefficients in (57) have reinforced symmetries which imply that

A+1,+1 = A-1,-1 and A+1,-1 = Al1,+1-

The right-hand side of Eq. (54) can be reduced further by expressing the spherical harmonics in terms of Legendre functions (Edmonds, 1960; Varshalovich et al., 1988). We have

i i Al Yiv-v(k b) Pi (cos 9 )

_ i 121 + 1

(- 2)!

+ . 7-rf 2 Re

(+ 2)!

IZ+il2 + IZ-il2

A+1,+1

Z-iZ+V2*1

A+1,-1 P<2)(cos 9)

where 9 and $ are the polar and azimuthal angles of the direction

kb, respectively, P>(cos 9) are Legendre polynomials and P^2)(cos 9) are associated Legendre functions of order 2. Finally, the DCS (54) can be cast in the form

d^WC) f i Ai Pi(cos 9 )

dkb W n6b 2IKa\ I ~ v 7

Re(Z-1ZÎ1) cos(2$) - lm(Z-iZ+i) sin(2$)

x I BiP<2)(cos 9)

where the coefficients A^ and B^ are real. It is worth noticing that

the DCS is invariant under inversion of the polarization vector

(Z ^ - Z). Evidently, expressions of the coefficients A^ and Bf are quite involved and hard to compute because of the abundance of vector coupling coefficients. In practice, it is more efficient to calculate the DCS from the less elaborate formula (49) and the expansion (50) of the emission amplitude.

k A I-

For linear polarization along the x axis, Z _ x ( Z+i _ - 1/V2, Z-i _ 1/V2 ), the DCS takes the form obtained by Scofield (1989):

d<n (W, x) (2n)2e2ch kb qa

W K6b 2|Ka|

x i { Pi(cos 9)- Bi cos(2$)P®(cos 9)}.

In the case of circular polarization, which corresponds to Z _ -( Z+1 _ - 1, Z-i _ 0, right handed) or Z _ ( Z+i _ 0, Z-i _ 1, left handed), Eq. (61) becomes

d<n(W, + l±i) (2n)2e2ch kb qa

W K6b 2\Ka\

i Ai Pi (cos 9).

As expected, in this case the angular distribution of photoelectrons is axially symmetrical, i.e., independent of the azimuthal angle <p

of the direction kb. We also note that the DCSs for right- and left-handed circular polarizations are identical, and equal to the DCS for unpolarized photons (see Section 3.1 ).

Sauter (1931) derived an analytical expression of the DCS for ionization of the K shell of atoms by photons linearly polarized along the x-axis from the Born approximation, i.e., from the above formulation with hydrogenic initial orbital and the photoelectron distorted plane waves replaced by Dirac plane waves. The normalized probability distribution function of the photoelectron direction obtained by Sauter reads (Sauter, 1931; Bethe and Salpeter, 1957)

Psauter (cos 9,

doK™ (W, x)

4on (W )

i + 3 Y (y - 2)

4 y + 1

sin2 9

(1 - fi cos 9)4

cos2 $

2fiY2 1 - fi

1 - 1Y (Y - 1)(1 - P cos 9 )

+ 1Y (Y - 1)2 (1 - P cos 9)

where fî and y are defined by Eqs. (32) and (33), respectively. The angular distribution for unpolarized photons is obtained by averaging over polarization directions (or, equivalently, over the azimuthal angle $), which gives

Psauter (cos 9)

%KY 4 sin2 9

i + 3 Y (Y - 2) 4 y + 1

2fiY2 1 - P

(1 - fi cos 9)4

1 + 2 Y (Y - 1)(y - 2)(i - fi cost

Sauter's distribution is used in the majority of high-energy Monte Carlo radiation codes. It gives fairly realistic results for the K shell of elements with small and moderate atomic numbers, say up to Z ~ 30. However, the corresponding total cross section is less accurate than Pratt's high-energy formula (31).

3.1. Partially polarized photons

Pure polarization states can be described by the Poincare vector, P = (P1, P2, P3), whose components are the Stokes parameters (see Appendix B). The Poincare vector of pure states has unit length and can be expressed in polar form [see Eq. (B.8)]

P = (sin a cos p, sin a sinp, cos a). (66)

The corresponding polarization vector is [see Eq. (B.10)]

Z(P) = cos(a/2)«Ai + sin(a/2) exp(ip)«A2, (67)

with spherical components [ = ZV Z, see Eqs. (A.35) and (A.33)]

•Se b m-27

_ Ar (Z= 18) - K shell

W (eV)

Fig. 2. Cross section for ionization of the K shell of argon atoms by absorption of photons, as a function of the photon energy W. The solid curve, which extends from threshold up to about 2 MeV, is the numerical result calculated from expression (16). The dashed curve (red online) is the prediction of the extrapolation formula (41). The dotted curve (blue online) was obtained from the dipole approximation. Diamonds represent data from the EPDL (Cullen et al., 1997).

Z+1(P) = - -L[cos(a/2) + i sin(a/2) exp(ifl)]],

Z-i(P) = -= [cos(a/2) - i sin(a/2) exp(ifl)]. It follows that

Z-i(p)Z*i(p) 2y—~ .......2

Consequently, the DCS (61) can be written as

d<" (W, P) (2^)2e2ch kb q,

- 1 (cos a - i sin a cos fl) = - 1 (P3 - iP1).

(68 a)

W neb 2|Ka|

z Al Pli

- [P3 cos(2^) + P1 sin(2^>)] Z BlPl2)(cos 9)

Let us now consider the angular distribution of photoelectrons for photons with arbitrary polarization represented by the Poin-care vector, P = (P1, P2, P3). The polarization density matrix can be expressed as [see Eq. (B.15)]

„ = (1 - P)1 i1 0J + P 1 i 1 + P3 P1 - iP2l, y y '2 {0 1J 2 {Pi + iP2 1 - P3 J

where P = (P12 + P| + P2)1/2 is the degree of polarization and P' = Pj/P. The latter quantities are the Stokes parameters of a pure state with polarization vector

Z(P') = cos(a'/2)1 + sin(a'/2) exp(ifl')$2, where

P1 + iP2

a ' = arccos P3,

exp (ifl ') :

41 - P3:

Eq. (71) means that the photon beam can be regarded as the incoherent superposition of two partial beams: an unpolarized beam with relative intensity (1 - P), and a fully polarized beam with

Poincaré vector P' = P/P and relative intensity equal to the degree of polarization. On the other hand, the unpolarized beam is equivalent to the superposition of two polarized beams with opposing polarizations, P' and -P', and equal intensities. This implies that the DCS for photons with polarization P can be calculated as a weighted average of the DCSs for photons with pure polarizations P' and -P',

dâ™ (W, P)

dk b 1 - P

dâ™ (W, Z (P ')) dâ™ (W, Z (- P ' ))

dkb dkb dâr (W, Z (P' ))

= dâa°n (w , z (p ' ))__

= o A +0 A

2 dkb 2 dkb

Z ( - P') = sin (a '/2) - cos(a '/2) exp (ifl')ê2

p dâr (W, Z (- P '))

is the polarization vector corresponding to the Poincare vector -P [see Eq. (B.10)]. The DCS (74) can be expressed as a Legendre series in terms of the Stokes parameters:

dâr (W, P)_(2n)2e2ch kb qa

W neb 2-al

Z Ai Pl(cos 9)

- [P3 cos(2^) + P1 sin(2^)] Z B lPl(2)(cos 9)

4. Numerical calculation of cross sections

Our Fortran program photacs calculates cross sections for ioni-zation and excitation of subshells of free atoms and ions using the theory presented above (Section 2). The structure and numerical

1 I 1 I 1 I 1 I 1 I 1 I 1

Ar (Z = 18) K shell

3175 3177 W(eV)

3177.4

3177.6 W (eV)

3177.8

Fig. 3. Left: cross section for absorption in the K shell of argon atoms as a function of the photon energy W in a narrow interval around the edge energy. The solid curve in the left panel is the cross section given by photacs-pp, including contributions from excitations to bound levels and accounting for the finite level width. Diamonds represent the energies and relative contributions (not to scale) of the lowest excited levels, which are treated as discrete resonances. Right: Magnified view of the step function (25), the polynomial approximation (29) (solid curve, red online) to the excitation quasi-continuum, and a part of the photoionization cross section (dashed curve). (For interpretation of the references to colour in this figure caption, the reader is referred to the web version of this paper.)

algorithms adopted in the program are similar to those employed in the programs developed by Bote and Salvat (2008) and Salvat and Bote (2011) to compute cross sections for inelastic collisions of charged particles with atoms. photacs calculates the radial functions of initial and final electron states by means of an updated version of the subroutine package radial (Salvat et al., 1995), which allows strict control of numerical round-off errors. The radial subroutines are also used in the self-consistent calculation of the DHFS potential. Vector-coupling coefficients and reduced matrix elements of the Racah tensors are calculated from their analytical formulas (Rose, 1995; Varshalovich et al., 1988) using a subroutine package that performs arithmetic operations at a high level of precision, well beyond Fortran double precision, by working in radix (base) 1000.

The program reads the potential V(r) felt by the active electron from an input file; the information to be provided in that file is a table of values of the function rV(r), for a dense enough grid of radii to allow accurate interpolation by natural cubic splines. The function rV(r) is required to be finite for all r, but is otherwise arbitrary. optionally, the program also allows the unscreened coulomb potential of the bare point nucleus to be used. As indicated above, in the calculations we normally use the DHFS potential described in Section A.2. For the evaluation of the integrals (A.42), the radial functions are calculated for a non-uniform radial grid, with 16 points in a wavelength, from which the integrals are evaluated using the 6-point Lagrange quadrature formula. The stability of the calculations was verified by using denser radial grid with 25 points/wavelength; cross sections computed with the two radial grids differ by less that 10-6.

The reduced cross section for excitation, agi0 (W), is computed according to Eq. (10). The difficulty of the calculation increases with the principal quantum number nb of the final level because the radial functions of the final state extend to larger radii. We only need values of the radial functions of the final level for computing the radial integrals (A.42) for the relatively small radii at which the radial functions of the initial state take appreciable values. However, the radial equations of the final state still need to be integrated up to large radii to determine the energy and normalization of the state. As indicated in Section 2.3, the program effectively computes the reduced cross sections for excitations to discrete levels with nb < 18.

The calculation of ionization cross section ol°n (W), Eq. (16), is performed by adding the grouped contributions from Kb = ± |Kb| in increasing order of \Kb|. The summation is discontinued at the first term which becomes less than 10-6 times the accumulated sum. Our program allows values of |Kb| to be considered up to ~200. The number of terms NK needed to get convergence of the series (16) increases with the energy of the photon, and it is fairly independent of atomic number of the target atom and the considered subshell. Typical values of N|K are about 5, 7,10,15, 35, and 150 for photons with energies of 100 eV, 1 keV, 10 keV, 100 keV, 1 MeV, and 10 MeV, respectively. An analysis of the variation of the calculated § ion(Kb, Ka) values with |Kb| indicates that the relative numerical errors of the computed ionization cross sections are usually less than about 10-5.

The ionization cross section &a°n (W) can be calculated for photon energies W from threshold (W = Ea) up to a cut-off energy WCut for which the spacing of the radial grid where wave functions are tabulated is insufficient to reproduce the fast oscillations of the integrands in Eqs. (A.42). Typically, the cut-off energy is larger than about 500Ea, where Ea = - ea is the ionization energy of the active shell. it is worth mentioning that, despite our superior computer power, photacs does not allow much higher energies than Scofield's calculations to be reached.

The computer program generates a table of &a°n (W) for a grid of energies extending from threshold up to the cut-off energy. This

bS 10-24

......... LI,2,3 Ar (Z= 18) -

- z2^ k K -

_ Ml,2,3 -

......... ......... .......i i ii 1"vsi'-j>o.i i

10J 104 W (eV)

Fig. 4. Cross sections for photoabsorption (excitation and ionization) in the subshells of argon atoms as functions of the photon energy W (thin solid curves), calculated by photacs-pp using experimental level widths. The thick solid curve (red online) represents the total atomic cross section. Dashed curves (blue online) are cross sections obtained from the dipole approximation. (For interpretation of the references to colour in this figure caption, the reader is referred to the web version of this paper.)

Fig. 5. Differential cross sections for ionization of the K shell of argon atoms by unpolarized photons of the indicated energies. Solid curves represent partial-wave numerical results. Dashed curves (blue online) are results from the dipole approximation. Dot-dashed curves (red online) correspond to the Sauter distribution, Eq. (65), multiplied by the numerical total cross section calculated from Eq. (16). (For interpretation of the references to colour in this figure caption, the reader is referred to the web version of this paper.)

table contains a nearly logarithmic grid, with 15 points per decade, plus a number of additional energies that are set by means of a self-adaptive method. The resulting energy grid is such that linear log-log interpolation errors are kept below some prescribed limit, which in the program was set equal to 0.05 %. For elements with atomic numbers near 19, 37, 55, and 87 (i.e., near the alkalies), this procedure reveals near-edge structures that were partially overlooked in Scofield's tables. otherwise, our results agree closely with those of Scofield, typically to within about 0.05 %, the differences being mostly attributable to the different numerical methods employed to solve the radial Dirac equations.

Fig. 2 displays the cross section for photoionization (W > Ea) of the K shell of argon atoms as a function of the photon energy. The solid line represents results calculated from the general formula (16), and extends over the energy interval where calculations with photacs are feasible, i.e., from the ionization threshold up to W^"1, which in the present case is 2 MeV. For comparison purposes, we

Fig. 6. DCSs for ionization of the K shell of Ar atoms by absorption of 50 keV photons, in barn (=10-24 cm2), as functions of the polar angle 9 and the azimuthal angle $ of the direction of the emitted photoelectron, both in degrees. The plots correspond to photons linearly polarized along the x-axis [ P = (0,0,1), top], along the y-axis [ P = (0, 0, - 1), middle], and along the direction Z = (x + y)/V2, at 45° from the x-axis [ P = (1, 0, 0), bottom].

have included cross-section data from the EPDL (Cullen et al., 1997), which were calculated by Scofield using the DHFS potential. The close agreement between our results and those of Scofield, which were generated using different computer codes, provides a clear indication of the accuracy of the numerical algorithms. Notice that the extrapolation formula (41) agrees well with the numerical results in a wide interval, differences are less than 2% down to about 200 keV.

Because excitation to bound levels is limited to a narrow energy interval, it is neglected in most of the existing databases and Monte Carlo simulation codes. Furthermore, the finite width of atomic levels is usually disregarded. These two effects become significant for describing the penetration and dosimetry of photon beams with energies near absorption edges. A situation where these effects are observed is in experimental measurements of the photon mass energy-absorption coefficient of air (Buermann et al., 2006). The values of the coefficient predicted by Monte Carlo codes, although in fairly good agreement with experiment, fail to reproduce the structure displayed by the measured coefficient near the energy of the K absorption edge of argon.

The contribution from excitation to bound levels (Section 2.3), the effect of the atomic-level width (Section 2.5), and the normalization screening correction (Section 2.6) are accounted for by a post-processing program, named photacs-pp, which reads the tables of numerical cross sections generated by photacs. This program determines the polynomial approximation (29) for the excitation pseudo-continua, extrapolates the ionization cross section to high energies using the formula (41), and performs a convolution with the Lorentzian profile (42) [see Eq. (44)]. The result from photacs-pp is a realistic photoelectric cross section, which varies continuously with energy and exhibits excitation structures that are in qualitative agreement with measurements in gases. Fig. 3 displays this cross section for absorption in the K shell of argon atoms of photons with energies near the edge (using the original DHFS energies calculated by photacs). The right panel is a magnified view of the step function (25) which describes the contribution of excitations to bound levels with nb < 18, and the polynomial approximation (29) to the quasi-continuum. Notice that, as we have already mentioned, the excitation and ionization cross sections match at the edge. Fig. 4 shows calculated cross sections for the subshells of argon and their sum, the atomic cross section. Interestingly, the dipole approximation is seen to predict subshell cross sections quite accurately for photon energies up to about 10 keV.

Optionally, the program photacs can calculate the DCS for ioni-zation, i.e., the photoelectron angular distribution, for a photon beam with energy W and polarization defined by the Stokes parameters P [see Eq. (74)]. The calculation starts by computing the ionization cross section ol°n (W), Eq. (15), as described above. The quantities e,m^<JtK6;JaKa [see Eqs. (53)] computed at this stage are stored in memory. The DCSs for pure polarization states can be calculated either from Eq. (49), which involves a minimum of angular momentum algebra, or from the more elaborate expression (61), which requires lengthier preliminary calculations. The

DCS is calculated and tabulated for a dense grid of directions kb of the photoelectron, defined by the polar and azimuthal angles, 9 and respectively. The DCSs obtained from the two schemes [i.e., from Eqs. (49) and (61)] are identical because we are using exact (double precision) values of vector coupling coefficients. The amount of work needed to sum the series (57) increases quickly with photon energy. Consequently, for photon energies well above the absorption edge, the calculation from Eq. (49) is much faster than from Eq. (61). The accuracy of the results can be verified by comparing the cross section ol°n (W) calculated previously with the value obtained by numerical integration of the DCS table for non-polarized photons. The numerical values resulting from the

two calculations normally agree to more than 5 digits.

Fig. 5 shows calculated DCS for ionization of the K shell of argon by unpolarized (or circularly polarized) photons. For comparison purposes, we have also included results from the dipole approximation, i.e., DCS calculated from the emission amplitude (50) with the G± coefficients given by Eq. (22). It is seen that the dipole approximation works much better for the total cross section than for the DCS. The dipole DCS agrees reasonably with the partial-wave results only for energies near the ionization threshold; it fails to describe the progressive decrease of the most probable emission angle when the photon energy increases. Interestingly, the Sauter distribution, Eq. (65), renormalized to reproduce the calculated cross section, provides quite an accurate description of the angular distribution for all energies because the example falls within the domain of applicability of Sauter's theory (K shell, moderate atomic number). For other subshells and for heavier elements the approximation is much less satisfactory.

Fig. 6 displays DCSs of the K shell of argon for 50 keV photons with linear polarizations along the x and y axes. The sum of these two distributions, with weights equal to 1, is the DCS for an unpolarized photon beam [see Eq. (74)], the numerical results agree with those obtained for a circularly polarized beam, which are independent of the azimuthal angle $ [cf. Eq. (63)], to more than 5 digits. The DCS for linear polarization at 45° from the x-axis, which corresponds to the Poincare vector P =(1,0, 0), is also displayed. Evidently, the maxima in the DCS are at directions with the azimuthal angle coinciding with that of the electric field; the polar angle of the maxima is close to 90° at small energies and decreases when the photon energy increases (Fig. 5) because photoelectrons do absorb part of the linear momentum of the photon.

5. Concluding comments

We have presented a detailed formulation of the theory of the photoeffect, within the one-active-electron approximation, in a form that is suitable for implementation in a computer program, complemented with extrapolation schemes to cover the full energy range of interest. Because of the robustness of the numerical methods (exact vector-coupling coefficients, highly accurate and densely tabulated radial functions), the program photacs provides reliable results for any atomic (or ionic) target and for any photon energy below the practical cut-offs. It allows generating cross section tables of subshell cross sections for the elements with unprecedented detail, including excitation to bound levels and angular distributions of photoelectrons for arbitrary photon polarizations. Although in the present calculations we have adopted the DHFS self-consistent potential, the program can also work with other screened potential models, e.g., to study the influence of the electron vacancy left after photoabsorption.

The post-processing program photacs-pp uses elaborate extrapolation models to determine subshell cross sections for highly-excited discrete levels and for photons with energies above the calculation cut-offs. It can also account for the effect of the finite width of atomic levels, which yields the cross section as a continuous function of the photon energy.

These programs have been used to generate a complete database of photoionization cross sections for the inner subshells (up to the N7 subshell) of the ground state configurations of the elements Z = 1-99, for photon energies from the ionization threshold up to 1 GeV. The provided cross-section tables were built from the numerical cross sections calculated by photacs (with the DHFS potential), which were extrapolated to energies higher than the calculation cut-off by means of Pratt's extrapolation formula,

shifted in energy to have the absorption edges coinciding with the experimental subshell ionization energies given by carlson (1975), and renormalized using McDF/DHFS density ratios. This database, which has already been adopted in the Monte Carlo code penelope (Salvat, 2015), is available under request.

Acknowledgments

We are indebted to A. Jablonski for useful discussions and advice on data required in X-ray photoelectron spectroscopy, and to J. E. Fernandez for kindly setting out our collaboration. Financial support from the Spanish Ministerio de Economía y Competitivi-dad and ERDF (Project no. FPA2013-44549-P) and from the Gen-eralitat de Catalunya (Grant 2014 SGR 846) is gratefully acknowledged. L.S. is thankful to the Department of Industrial Engineering of Alma Mater Studiorum University of Bologna for financial support through the Marco Polo Grant programme.

Appendix A. Wave functions and matrix elements

We describe here various notations and mathematical formulas used in the present calculations of the photoeffect. These calculations are based on a model of independent electrons in a central potential. in other words, states of the target atom are represented as single Slater determinants, built with orbitals that are eigen-functions of the one-electron Dirac Hamiltonian for a central potential V(r) (Rose, 1961)

"H = c ap + ß mec2 + V ( r),

where p = - ihV is the momentum operator, and a, ft are the Dirac 4 x 4 matrices in the spinor representation:

¡2 0 0 - ¡2

( a a), ¡>

The 2 x 2 matrices a are the Pauli matrices:

=( îi) =(0 V) =(1 -1}

A.1. Electron wave functions

The spherical orbitals are simultaneous eigenfunctions of the Dirac Hamiltonian and the total angular momentum J = L + S, where L = - ir x V is the orbital angular momentum and S is the spin angular momentum (all angular momenta in units of h). These eigenfunctions are the so-called spherical waves, and have the form (Rose, 1961; Grant, 1965)

Vexm ( r) = -

PeK ( r ) QK,m () iQex ( r)Q-k,m ()

where QKm (^ are spherical spinors, and PeK (r) and QeK (r) are the large- and small-component radial functions, respectively, which satisfy the coupled differential equations

dPjK _ k P e - V + 2meC2 Q

— = peK + 1 ' ~ -

' ex 1 a.

e-V P + k 0

where e is the electron energy, exclusive of its rest energy. The spherical spinors are eigenfunctions of the total angular

momentum of pauli's theory, i.e., simultaneous eigenfunctions of the operators L2, Sp, / and JPz with eigenvalues 1(1+ 1), 3/4, j (j + 1) and m, respectively. Here SP = i<7 denotes the two-dimensional Pauli spin operator and Jp = L + SP. We have (Rose, 1995; Varshalovich et al., 1988)

,«=±1/2

I, -1, m - «, « j, m

V,m- «

The quantities (jjm^m2 j, m) are Clebsch-Gordan coefficients, Ylm (r) are spherical harmonics, and ^ are the Pauli unit spinors, i.e., the eigenstates of Sp and SP3 with eigenvalues 3/4 and H = ± 1, respectively. More explicitly,

± ± m + I Yi,m-l/2 (r)

Q/±i/2m (r) =

11 + m + - Yi,m+1/2

To simplify notation, it is customary to use the relativistic angular momentum quantum number

K = (- j)(2j + 1)

which specifies both the total angular momentum, j, and the parity, ( - 1)1 of the Dirac spherical wave:

j = K - 2,

I = j + t;— 2 |K

K if K >

- K - 1 if K <

It is also convenient to consider the quantum number

- K if K < K - 1 if K >

= I - ■

(A.10)

which is the value of I corresponding to -k.

In the case of bound orbitals (e < 0), each discrete energy level is characterized by the principal quantum number n and the re-lativistic quantum number k. Bound orbitals calculated by the radial subroutines are normalized to unity and, consequently, satisfy the orthonormality relation:

y^K'm'( r) Wnxm ( r) dr = S n'n Sx k S „

(A.11)

Free spherical waves (with e > 0) are normalized in such a way that the large-component radial function asymptotically oscillates with unit amplitude

PeK(r) ~ sin( kr - fn - n ln 2kr + S I, V 2 J

k = (ch)-1Je(e + 2mec2)

(A.12)

(A.13)

is the wave number, n = rV(r)] me/(Ä2k) is the Sommerfeld

parameter, and öK is the phase shift. Free spherical waves normalized in the form (A.12) satisfy the orthogonality relation

Wex-m'( r) Vex

,( r) dr = ^ n S (e' - e) Sx k S mm

(A.14)

The state of a free electron with definite spin projection can be represented as a distorted plane wave (DpW), i.e., a solution of the Dirac equation that asymptotically behaves as a plane wave plus an outgoing (+) or incoming (-) spherical wave. A DPW is characterized by the wave vector k of the asymptotic plane wave and the spin projection ¡i; it can be expanded in the basis of spherical waves as (see, e.g., Rose, 1961; Pratt et al., 1973; Akhiezer and Berestetskii, 1965)

<( r>= 1

^ kv n(e + mec2

£ il exp( ± i5K) [QKm(k)]

X VeKm ( r) •

e = J(chk)2 + (mec2)2 - mec2

(A.15)

(A.16)

is the kinetic energy of the particle. The expansion (A.15) is known as the partial-wave series. It can be easily verified that, with the adopted normalization for free spherical waves, DPWs satisfy the orthogonality relation

/[ vk±kr)]Vk±(r) dr = 8(k' - k) (A.17)

A.2. The DHFS potential

The DHFS potential of an atom or ion of atomic number Z with N bound electrons is determined by the ground-state electron density, p(r), which is obtained self-consistently (see, e.g., Liber-man et al., 1965, 1971). It is given by

VdHFS (r) = V"uc (r) + Vei ( r) + Vex ( r).

(A.18)

The term Vnuc (r) is the nuclear potential (i.e., the electrostatic interaction energy of an electron at r with the nucleus, which is assumed to be spherical). For a point nucleus,

Vnuc(r)= - Ze2/r.

(A.19)

The effect of the finite size of the nucleus can be accounted for by using simple models for the nuclear charge distribution. A convenient parameterization of the proton density is provided by the Fermi distribution (Hahn et al., 1956)

Pp(r ) =

exp[(r - Rn)/z] + 1'

Rn = 1.07 Zl1/3fm, and z = 0.546 fm,

(A.20)

(A.21)

where A is the mass number, which is usually replaced by the atomic weight (mean relative atomic mass) of the element. The constant p0, which equals twice the proton density at r = Rn, is to be determined by normalization. The nuclear potential for the Fermi distribution has to be calculated numerically:

Vnuc(r

e2 [fp^l dr' J |r - r 'l

e2 pv /»to

-7 /0 Pp(r,) 4nr'2 dr' - e2 /r PP(r')4nr' dr' (A.22)

The second term in expression (A.18) is the electronic potential (i.e., the interaction energy of an electron at r with the atomic electron cloud)

e2 pr /»to

Vei(r) = — / p(r')4nr'2dr' + e2 / p(r')4nr' dr', r J0 Jr

(A.23)

and the last term is the local approximation to the exchange potential. We use the Slater-Latter potential given by

VeSXater (r)

if r < rLatter,

Vex (r )= < e2

(Z - N + 1)— - Vnuc (r)- Vei (r) if r > rLatter, (A 24)

VeS,xater (r) = - 3e2 (3/n)1/3 [p (r)]1/3.

(A.25)

is the exchange potential derived by Slater (1951), and the cut-off radius rLatter is the outer root of the equation:

r [ Vnuc(r) + Vei(r) + aV|xater(r)] = - (Z - N + 1)e2

(A.26)

The modification (A.24) of Slater's potential for r > rLatter, which ensures the correct behaviour of V(r) at large radii, is known as Latter's tail correction (Latter, 1955).

A.3. Matrix elements ofRacah tensors

The calculation of the matrix elements (3) and (7) is performed by expanding the plane wave in terms of the Racah functions

CLM (A) - j^TT YLM

(A.27)

to allow analytical evaluation of the spin and angular parts. We recall that the 2L + 1 functions CLM(r) constitute an irreducible tensor of rank L, C(L). By virtue of the Wigner-Eckart theorem (see, e.g., Edmonds, 1960), the matrix elements of Racah tensors for eigenstates QKm of the total angular momentum J = L + SP of a spin 1 particle are

{QK1m1 |Qm 2m2) - / [ ^K1m1(()]|t Clm ( () 2m2 (() dr

= 72171 (j2 Lm2Mlj1m1)

x(f1j\ C^^). (A.28)

The reduced matrix element ^li-jji ||C(L)||l21j2 ^ is given by the expression (see, e.g., Grant, 1961)

C(L)K->2 ) = U (L, 11, l) j+T (LJ201| j1jj,

where the factor

v (L, 11, 12)

1 if L + 11 + l2 is even, 0 otherwise,

(A.29)

(A.30)

accounts for the parity selection rule.

Matrix elements of spin operators can be evaluated by using the identity

^Kt1m1 (r)ff 2m2 (

= £ (-1)M u j

j k + k2

1K ^^k 1,m1 |Cjm |fl-K 2,m2) YJ--M (A)

j^K +1 1) (k1 - K2) (fl^m,|CJMpK2,m2) YJ,-M(A)

jr ( 7+T + 1}Q*m\CmP-™) jm(A)

(A.31)

where the functions YLJM(r) are the vector spherical harmonics, defined by (see, e.g., Varshalovich et al., 1988; Akhiezer and Be-restetskii, 1965)

YLM(() = £ (l, 1, M + v, - v| JM) Ylm+v(() i-

(A.32)

where IV are the spherical unit vectors:

(A.33)

Calculations are easier when vectors G are expressed in the spherical basis

G = X ( - 1)"^Ç-v

where Gv _ Çv G

(A.34)

(A.35)

are the spherical components, which constitute an irreducible tensor operator (Rose, 1995; Edmonds, 1960). The dot product of two vectors, Z and G, is given by

Z G _ x ZiGi _ x( - 1)vZvG-v,

(A.36)

where Q, Gj and Gv are, respectively, the Cartesian and spherical components of the vectors.

A.4. One-electron transition-matrix elements

The matrix elements (3) and (7) are of the type

Mba = Z G with G = (yebKbmb |a exp (ik-r)^eaKam^ , (A.37)

where, in general, the polarization vector Z may be complex. To calculate the spherical components of G, Eq. (A.35), we introduce the Rayleigh expansion of the plane wave

to j . a.

exp (ik-r) = x x iJ(2J + 1))J(kr)CM(k)(if),

J=0 M=-J

(A.38)

where jj (kr) are spherical Bessel functions. The formulas become simpler if we consider a reference frame with the z-axis parallel to the direction k of the photon. In such a frame, CfM (k) = ôM0, and we can write

= x ij (2j + 1)(webkbmb\a jj (kr)Cjo (r)LeaKama\.

(A.39)

The matrix elements in this expression can be evaluated on the basis of Eq. (A.31), by a method similar to the one adopted by Mann and Johnson (1971) for a related problem. After a rather tedious calculation we find the following expressions for the spherical components of G

G+1((, mb; Ka, ma )

to r -I

= x iJ 2J + 1

pJ ( + 1)

{ (^Kb,mb|Cj,+1 |^Ka,ma) e^ebKb;eaKa

i (^Kb,mb|CJ,+1 |^-Ka,ma) ^ebKb;eaKa},

(A.40a)

Ebk b; Eaka

with the radial integrals

= 2J + 1 [ (( - Ka) (bKb;eaka + Gebkb;eaka )

- J (FebXb;e + (( - Ka) (( ;e + (( + 1)(^!j+Kb;eaKa - Gebkb;eaka

- GJ-1 )

,;eaKa ^ebk-b;eak-a)

Gebkb;eaka)

(A.41a)

J( + 1

a»J _ _

ebk-b; eaka 2J + 1

Kb - Ka

fj-1 + GJ-1

ebkb;eaka ^ uebkb;eaka

FJ-1 - Gj-1

rebkb;eaka Kjebkb;eaka

Kb - Ka I Fj+1 + Gj+1

I ebkb;eaka T uebKb;eaKa

- I FJ+1 - GJ+1

1 1ebkb;eaka uebkb;eaka

<3bKb;eaKa _ ( Ka + Kb) (F4Kb;eaKa + Glbk-b;eak-a ), where

/•to

FeWaKa = PebKb (r) QeaKa (r) jj (kr)

ebK-b;eaK a

( r ) PeaKa (r) jj (kr) dr.

(A.41b)

(A.41c)

(A.42a)

(A.42b)

The superscripts "l", "e" and "m" in the integrals (A.41) stand for "longitudinal", "electric" and "magnetic", respectively, because these radial integrals also arise in an alternative treatment based on the multipole expansion of the radiation field (see, e.g., Scofield, 1975, 1978).

We recall that the matrix elements (QKb,mb\Cj,vQKa,ma) vanish unless mb = ma + v. Moreover, at least one of the elements (^K6,mb|C;,v|^Ka,ma) and Q^Cj,v|^-Ka,ma) vanishes because the values of la for Ka and -Ka differ by one unit and, therefore, the parity factor u (J, lb, la) is null for one of these matrix elements. Hence, the two terms in curly braces on the right-hand side of Eqs. (A.40a) and (A.40b) cannot be both different from zero.

More compact formulas are obtained by introducing the quantities

' Jb;eK - JT {fb j || C(L)|la j )

ebkb;eaka ■

2J + 1

T)(l j IMK1'a}

(A.43a)

(A.43b)

G-1(b, mb, Ka, ma) to ") J 1 _ x ij { (^Kb,mb|CJ,-1 |^Ka,m^ <

- i (^Kb,mb|CJ,-1 |^-Ka,m^ m^ebKb;eaKa},

Go((, mb; Ka, ma)

_ x iJ (2J + 1H«Kb ,mb Cj,0 ^^Ka ,mal ^ebKb;eaKa<

(A.40b)

(A.40c)

ebkb;eaka '

Jj "+11)lb IMI4 :Tja } m^JJ^b;c

(A.43c)

where la is the orbital angular momentum quantum number corresponding to -Ka. Considering the symmetry properties of Clebsch-Gordan coefficients we can write

G±1 _ ( - 1)ja-ma X iJ j, jb, ma, - mb |j, + 1)

± i myJ

ebkb;eaka - u ebkb;eaka

(A.44a)

G0 = ( - 1)ja-ma £ iJj, jb, ma, - mb\J, 0} XJ J=0

(A.44b)

The cross section for photoabsorption by the electrons in a closed subshell (naKa)2lKal is proportional to the sum of squared transition-matrix elements over magnetic quantum numbers, Eqs. (10) and (17),

gba = £ Mba P ma,mb

Cv> £ Gv G *

(A.45)

Since the photon polarization vector Z is perpendicular to the wave vector k = kz and normalized to unity, we have f0 = 0 and IC_i|2 + If+il2 = 1. We note that the only dependence on the quantum numbers ma and mb is through the clebsch-Gordan coefficients (ja, jb, ma, - mbJ, v). The summation over these quantum numbers is performed easily by using the orthogonality property of these coefficients. We thus obtain,

£ G±1G*1 = 0,

£ G±1G±1 = £ j [ ^^é6Kb;eaKa] + [ ^ebKb;eaKa] f ma,mb J=1 [

Therefore, the quantity (A.45) is given by

Qba = £ jf[e^^ébKb;eaKa] + [^^CaKa] j-

(A.46a)

(A.46b)

(A.47)

Appendix B. Photon polarization

The polarization of photon beams that propagate in the direction of the z-axis can be described in terms of the basis { \e2Sj } of linear polarization states along the x and y axes. Pure polarization states can be expressed in the form

|Z) = cos(a/2)(Z^ + sin(a/2) exp(i cos(a/2) sin(a/2) exp(ift) 0

with a e [0, n] and ft e ( - n, n]. There is a one-to-one correspondence between polarization states |Z) and polarization vectors

of electromagnetic waves, Z = S^i + Z2e2. The electric field £ of the polarized wave is given by

£ (r, t) = Re { [ cos (a/2)<A1 + sin (a/2) exp (ift )e2 ]

x exp[i(kr - mt + j)]}. (B 2)

In the following, we simplify the notation by writing Z and €[ instead of |Z) and , i.e., we use the same symbol to designate the quantum polarization states and the unit polarization vectors. States with ft = 0 correspond to linear polarization in a direction that makes an angle a/2 with e1. If ft = ± n/2 and a = n/2, we have right-handed (r) and left-handed (l) circularly polarized photons,1

1 We consider "natural" right-handed polarization, which is opposite to the convention adopted in optics (Born and Wolf, 2002).

e(r) e(D

= 7T( + *2)= - *+ = 12- i^2 ) = 1-1,

where and are the vectors of the spherical basis, Eq. (A.33).

Usually, real beams are partially polarized, that is, their photons are in various states |Zn ) of pure polarization with corresponding probabilities pn. The states |Zn ) are only assumed to be normalized to unity, their number and nature are arbitrary. The probabilities pn are positive and add to unity. polarization features of real beams can be described by using the density matrix formalism (Falkoff and MacDonald, 1951; Fano, 1954, 1957; McMaster, 1954). The density operator for such a radiation beam is

P = £ |Z„)Pn(Zn|-

The matrix of this operator in the basis of states of linear polarization {e1, e2} is Hermitian and has unit trace. It can be expressed as a linear combination of the Pauli matrices a, Eq. (A.3), and the 2 x 2 identity matrix, J2, with real coefficients:

p = 1 (/2 + + P2ff2 + P3*3) = 1 f 1 + P3 P1 - iP2l F 2yJ 1 1 2 2 3 3 2 {P1 + iP2 1 - P3 J (B.5)

The coefficients Pi are the Stokes parameters; they provide a complete description of the polarization of a beam and can be measured experimentally (Fano, 1954). We have Pi = Tr po\ or, more explicitly,

P1 = P12 + P21 P2 = i (P12 - P21), P3 = P11 - P22*

It is worth mentioning that in optics the pauli matrices, as well as the Stokes parameters, are usually considered in a different order, namely, {03, a1, o2}. We prefer the ordering (A.3) employed in quantum mechanics because the formalism is thus parallel to that of polarization of spin-1 particles.

The Stokes parameters can be regarded as the components of a vector, P = (P1, P2, P3), the Poincaré vector, which is analogous to the direction of spin of a spin-1 particle, although it transforms differently under rotations (see Fano, 1957, 1954). For a pure state of the type (B.1 ) the density matrix takes the form

2 (a/2) cos (a/2) sin (a/2) exp ( - ifl)

sin2 (a/2)

cos(a/2) sin(a/2) exp(ifl)

and the associated Stokes parameters are P1 = sin a cos ft, P2 = sin a sin ft, P3 = cos a.

Note that a and ft are the polar and azimuthal angles of the poincaré vector P, respectively, and that P= 1. In the opposite situation, when the Stokes parameters vanish, P=0, the density matrix takes the form P = 1/2, which represents unpolarized photons. The magnitude P of the poincaré vector measures the degree of polarization; it can take values from 0 (unpolarized photons) to 1 (pure polarization states).

In the case of pure states (P = 1), inverting the relations (B.8), we obtain the state angles (a, ft) from the Stokes parameters:

a = arccos P3, exp (ift) :

P1 + iP2

VT-Pf '

The pure states corresponding to the poincaré vectors P and with respective directions (a, ft) and (n - a, ft + n), are

(B.9) P,

cos(a/2) Z (P)= sin(a/2) exp (ift) and 0

sin(a/2) Z(- P)= - cos(a/2) exp(ift) 0

Note that these states are orthogonal:

(Z (P) IZ (- P) = 0.

1(1 + P3 Pi - 1P2

(B.10)

(B.11)

Hence, by reversing the signs of the Stokes parameters of a pure state, we obtain its orthogonal state (except for, possibly, an irrelevant phase factor). Thus, the state angles (a, ft), the Poincare vectors and the density matrices of the states {«A1, e2} of the linear-polarization basis are

C2 = Y =

a = 0 ft = 0'

a = n ft = n '

0 0 - 1

1 f 2 rt 2 f 0 0/

( 0 0).

(B.12a)

(B.12b)

Similarly, for the states of the basis of circular polarization, we have

l+1 = J2

_ a = n/2 ' ft = n/2'

if 1 v),

(B.13a)

1-1 = V2

a = n/2 ' ft = - n/2'

0 -1 0

1 f -1 1}

(B.13b)

The intensity of a photon beam is defined as the number of photons that cross a unit surface perpendicular to the direction of the beam per unit time. Let us consider two photon beams with intensities I1 and I2 and respective density matrices p1 and P2. The incoherent admixture of these two beams gives a beam with intensity / = /1 + /2 and density matrix

/1 ¡2

P = / P1 + j P2.

The corresponding Stokes parameters are

(B.14a)

(B.14b)

where P1 and P2 are the poincare vectors of the initial beams. Notice that an unpolarized beam can be considered as the admixture of two completely polarized beams with equal intensities and "opposite" polarizations, P and -P, whatever the direction of P.

A partially polarized beam with Stokes parameters P (P < 1) can be regarded as an incoherent admixture of an unpolarized beam and a completely polarized beam. To characterize these beams, we define the reduced Stokes parameters P' - p/P. The matrix density of the partially polarized beam P can then be expressed as

2 f P1 + iP2 1 - P3

=(1 - 11)+Pi

1 + P3 p1 + iP'

p1 - iP' 1 - P3

(B.15)

Hence, the original beam is equivalent to the mixture of an un-polarized beam and a completely polarized beam (with poincare vector P' = P/P, P' = 1), having relative intensities (1 - P) and P, respectively.

References

Akhiezer, A.I., Berestetskii, V.B., 1965. Quantum Electrodynamics. Interscience, New York.

Amusia, M.Ya., Baltenkov, A.S., Chernysheva, L.V., Felfli, Z., Msezane, A.Z., 2001. Nondipole parameters in angular distributions of electrons in photoionization of noble-gas atoms. Phys. Rev. A 63, 052506.

Baym, G., 1974. Lectures in Quantum Mechanics. Westview Press, Boulder, Colorado.

Berger, M.J., Hubbell, J.H., Seltzer, S.M., Chang, J., Coursey, J.S., Sukumar, R., Zucker, D.S., 2005. XCOM: Photon Cross Sections Database. Technical report, Institute of Standards and Technology, Gaithersburg, MD. Available from (http://physics. nist.gov/PhysRefData/Xcom/Text/XCOM.html>.

Bethe, H.A., Salpeter, E.E., 1957. Quantum Mechanics of One- and Two-Electron Atoms. Springer-Verlag, Berlin.

Born, M., Wolf, E., 2002. Principles of Optics, 7th (expanded) edition. Cambridge University Press, Cambridge.

Bote, D., Salvat, F., 2008. Calculations of inner-shell ionization by electron impact with the distorted-wave and plane-wave Born approximations. Phys. Rev. A 77, 042701.

Büermann, L., Grosswendt, B., Kramer, H.-M., Selbach, H.-J., Gerlach, M., Hoffmann, M., Krumrey, M., 2006. Measurement of the X-ray mass energy-absorption coefficient of air using 3 keV to 10 keV synchrotron radiation. Phys. Med. Biol. 51, 167-175.

Burke, P.G., Taylor, K.T., 1975. R-matrix theory of photoionization. Application to neon and argon. J. Phys. B: Atom. Mol. Phys. 8, 2620-2639.

Campbell, J.L., Papp, T., 2001. Widths of the atomic K-N7 levels. At. Data Nucl. Data Tables 77, 1-56.

Carlson, T.A., 1975. Photoelectron and Auger Spectroscopy. Plenum Press, New York.

Chang, T.N., Fano, U., 1976. Many-body theory of atomic transitions. Phys. Rev. A 12, 263-281.

Cooper, J.W., 1990. Multipole corrections to the angular distribution of photoelec-trons at low energies. Phys. Rev. A 42, 6942-6945, Erratum, Phys. Rev. A, 45 (1992) 3362.

Cooper, J.W., 1993. Photoelectron-angular-distribution parameters for rare-gas subshells. Phys. Rev. A 47, 1841-1851.

Cullen, D.E., Hubbell, J.H., Kissel, L., 1997. EPDL97 The Evaluated Data Library, '97 Version. Technical Report UCRL-50400, Lawrence Livermore National Laboratory, Livermore, California.

Derevianko, A., Hemmers, O., Oblad, S., Glans, P., Wang, H., Whitfield, S.B., Wehlitz, R., Sellin, I.A., Johnson, W.R., Lindle, D.W., 2000. Electric-octupole and pure-electric-quadrupole effects in soft-X-ray photoemission. Phys. Rev. Lett. 84, 2116-2119.

Desclaux, J.P., 1975. A multiconfiguration relativistic Dirac-Fock program. Comput. Phys. Commun. 9, 31 -45.

Desclaux, J.P., 1977. Erratum notice. Comput. Phys. Commun. 13, 71.

Edmonds, A.R., 1960. Angular Momentum in Quantum Mechanics. Princeton University Press, Princeton, NJ.

Falkoff, D.L., MacDonald, J.E., 1951. On the Stokes parameters for polarized radiation. J. Opt. Soc. Am. 41, 861-862.

Fano, U., 1954. A Stokes-parameter technique for the treatment of polarization in quantum mechanics. Phys. Rev. 93, 121-123.

Fano, U., 1957. Description of states in quantum meachanics by density matrix and operator techniques. Rev. Mod. Phys. 29, 74-93.

Fano, U., Cooper, J.W., 1968. Spectral distributions of atomic oscillator strengths. Rev. Mod. Phys. 40, 441-507.

Fernández, J.E., Hubbell, J.H., Hanson, A.L., Spencer, L.V., 1993. Polarization effects on multiple scattering gamma transport. Radiat. Phys. Chem. 41, 579-630.

Fujikawa, T., Suzuki, R., Arai, H., Shinotsuka, H., Kovér, L., 2007. Nondipole effects in photoemission angular distribution excited by high-energy x-rays. J. Electron Spectrosc. Relat. Phenom. 159, 14-23.

Gavrila, M., 1959. Relativistic K-shell photoeffect. Phys. Rev. 113, 514-526.

Grant, I.P., 1961. Relativistic self-consistent fields. Proc. R. Soc. A 262, 555-576.

Grant, I.P., 1965. Relativistic self-consistent fields. Proc. Phys. Soc. A 86, 523-527.

Hahn, B., Ravenhall, D.G., Hofstadter, R., 1956. High-energy electron scattering. Phys. Rev. 101, 1131-1142.

Hemmers, O., Guillemin, R., Lindle, D.W., 2004. Nondipole effects in soft X-ray photoemission. Radiat. Phys. Chem. 70, 123-147.

Hubbell, J.H., Gimm, H.A., 0verb0, I., 1980. Pair, triplet, and total cross sections (and mass attenuation coefficients) for 1 MeV-100 GeV photons in elements Z — 1 to

100. J. Phys. Chem. Ref. Data 9, 1023-1147. Johnson, W.R., Lin, c.D., 1979. Multichannel relativistic random-phase approximation for the photoionization of atoms. phys. Rev. A 20, 964-977.

Kover, L., 2010. X-ray photoelectron spectroscopy using hard x-rays. J. Electron Spectrosc. Relat. Phenom. 178-179, 241-257.

Latter, R., 1955. Atomic energy levels for the Thomas-Fermi and Thomas-Fermi-Dirac potential. Phys. Rev. 99, 510-519.

Lee, P.A., Citrin, P.H., Eisenberger, P., Kincaid, B.M., 1981. Extended X-ray absorption fine-structure—its strengths and limitations as a structural tool. Rev. Mod. Phys. 53, 769-806.

Liberman, D., Zangwill, A., 1984. A relativistic program for optical response in atoms using a time-dependent local density approximation. Comput. Phys. Commun. 32, 75-82.

Liberman, D.A., Cromer, D.T., Waber, J.T., 1965. Self-Consistent-field Dirac wave functions for atoms and ions. I. Comparison with previous calculations. Phys. Rev. 137, A27-A34.

Liberman, D., Cromer, D.T., Waber, J.T., 1971. Relativistic self-consistent field program for atoms and ions. Comput. Phys. Commun. 2, 107-113.

Mann, J.B., Johnson, W.R., 1971. Breit interaction in multielectron atoms. Phys. Rev. A 4, 41-51.

McMaster, W.H., 1954. Polarization and the Stokes parameters. Am. J. Phys. 22, 351-362.

Parpia, F.A., Johnson, W.R., 1983. Application of the time-dependent local density approximation to the photoionisation of mercury. J. Phys. B: Atom. Mol. Phys. 16, L375-L379.

Perkins, S.T., Cullen, D.E., Chen, M.H., Hubbell, J.H., Rathkopf, J., Scofield, J., 1991. Tables and Graphs of Atomic Subshell and Relaxation Data Derived from the LLNL Evaluated Atomic Data Library (EADL), Z = 1-100. Technical Report UCRL-iD-50400, Lawrence Livermore National Laboratory, Livermore, California.

Pratt, R.H., 1960a. Atomic photoelectric effect at high energies. Phys. Rev. 117, 1017-1028.

Pratt, R.H., 1960b. Photoeffect from the L shell. Phys. Rev. 119, 1619-1626.

Pratt, R.H., Tseng, H.K., 1972. Behavior of electron wave functions near the atomic nucleus and normalization screening theory in the atomic photoeffect. Phys. Rev. A 5, 1063-1072.

Pratt, R.H., Ron, A., Tseng, H.K., 1973. Atomic photoelectric effect above 10 keV. Rev. Mod. Phys. 45, 273-325.

Richtmyer, F.K., Barnes, S.W., Ramberg, E., 1934. The widths of the L-series lines and of the energy levels of Au(79). Phys. Rev. 15, 843-860.

Rose, M.E., 1961. Relativistic Electron Theory. John Wiley and Sons, New York.

Rose, M.E., 1995. Elementary Theory of Angular Momentum. Dover, New York.

Sakurai, J.J., 1967. Advanced Quantum Mechanics. Addison and Wesley, New York.

Saloman, E.B., Hubbell, J.H., Scofield, J.H., 1988. X-ray attenuation cross sections for energies 100 eV to 100 keV and elements Z= 1 to Z= 92. At. Data Nucl. Data Tables 38, 1-197.

Salvat, F., 2015. penelope-2014: A Code System for Monte Carlo Simulation of Electron and Photon Transport. OECD/NEA Data Bank, NEA/NSC/D0C(2015)3, Issy-les-

Moulineaux, France, 2015. Available from <http://www.nea.fr/lists/penelope. html).

Salvat, F., Bote, D., 2011. Inelastic Collisions of Fast Charged Particles with Atoms. Relativistic Plane-wave Born Approximation. Technical Report, Universitat de Barcelona, Spain.

Salvat, F., Fernández-Varea, J.M., 2009. Overview of physical interaction models for photon and electron transport used in Monte Carlo codes. Metrologia 46, S112—S138.

Salvat, F., Fernández-Varea, J.M., Williamson, W., 1995. Accurate numerical solution of the radial Schrodinger and Dirac wave equations. Comput. Phys. Commun. 90, 151-168.

Sauter, F., 1931. Über den atomaren Photoeffekt in der K-Schale nach der relativistischen Wellenmechanik Diracs. Ann. Phys. 11, 454-488.

Schmickley, R.D., Pratt, R.H., 1967. X-, L-, and M-shell atomic photoeffect for screened-potential models. Phys. Rev. 164, 104-116.

Scofield, J.H., 1973. Theoretical Photoionization Cross sections From 1 to 1500 keV. Technical Report UCRL-51326, Lawrence Livermore Laboratory, Livermore, California.

Scofield, J.H., 1975. Radiative transitions. In: Crasemann, B. (Ed.), Atomic Inner-Shell Ionization, Ionization and Transition Probabilities vol. 1. Academic Press, New York, pp. 265-289.

Scofield, J.H., 1978. K- and L-shell ionization of atoms by relativistic electrons. Phys. Rev. A 18, 963-970.

Scofield, J.H., 1989. Angular and polarization correlations in photoionization and radiative recombination. Phys. Rev. A 40, 3054-3060.

Slater, J.C., 1951. A simplification of the Hartree-Fock method. Phys. Rev. 81, 385-390.

Trzhaskovskaya, M.B., Nefedov, V.I., Yarzhemsky, V.G., 2001. Photoelectron angular distribution parameters for elements Z— 1 to Z=54 in the photoelectron energy range 100-5000 eV. At. Data Nucl. Data Tables 77, 97-159.

Trzhaskovskaya, M.B., Nefedov, V.I., Yarzhemsky, V.G., 2002. Photoelectron angular distribution parameters for elements Z—55 to Z —100 in the photoelectron energy range 100-5000 eV. At. Data Nucl. Data Tables 82, 257-311.

Trzhaskovskaya, M.B., Nikulin, V.K., Nefedov, V.I., Yarzhemsky, V.G., 2006. Non-di-pole second order parameters of the photoelectron angular distribution for elements Z —1-100 in the photoelectron energy range 1-10 keV. At. Data Nucl. Data Tables 92, 245-304.

Varshalovich, D.A., Moskalev, A.N., Khersonskii, V.K., 1988. Quantum Theory of Angular Momentum. World Scientific, Singapore.

Watanabe, T., 1965. Measurement of the L absorption spectra of xenon. Phys. Rev. 137, A1380-A1382.

Williams, G.P., 2011. Electron binding energies of the elements. In: Haynes, W.M., Lide, D.R. (Eds.), CRC Handbook of Chemistry and Physics, 91st edition CRC, Boca Raton, pp. 221 -226.

Zangwill, A., Soven, P., 1980. Density-functional approach to local-field effects in finite systems: photoabsorption in the rare gases. Phys. Rev. A 21, 1561-1572.