A theoretical model of delayed emission following photoexcitation from metals and semiconductors is given. Its numerical implementation is designed for beam optics codes used to model photocathodes in rf photoinjectors. The model extends the Moments approach for predicting photocurrent and mean transverse energy as moments of an emitted electron distribution by incorporating time of flight and scattering events that result in emission delay on a subpicosecond level. The model accounts for a dynamic surface extraction field and changes in the energy distribution and time of emission as a consequence of the laser penetration depth and multiple scattering events during transport. Usage in the ParticleinCell code MICHELLE to predict the bunch shape and duration with or without laser jitter is given. The consequences of delayed emission effects for ultrashort pulses are discussed.
I. INTRODUCTION
Beam optics codes are essential for the successful prediction of quantum efficiency QE and emittance $ \epsilon n , rms $ from photocathodes and for the design of particle accelerators, next generation light sources, Free Electron Lasers (FELs), and other applications requiring emittance dominated beams.^{1–4} When the laser pulse is long compared to the response time of the photocathode, the reliance on instantaneous emission models is acceptable and widely done, but for future xray Free Electron Lasers (xFELs) which require pulse lengths to be compressed to, or when pulse shaping demands, subpicosecond risefall times, the approximation is no longer accurate. This is particularly true for semiconductors that have deep laser penetration effects and for which phonon scattering mechanisms do not reduce photoexcited electron energy as quickly as for metals: both contribute to emission tails.^{5–8}
The Moments model of photoemission^{9–11} (and the closely related Spicer ThreeStep model^{12,13}) treats photoemission by considering three processes: (i) absorption of the photon by electrons that have accessible final states; (ii) transport by the electron to the surface with losses due to scattering; and (iii) emission into vacuum by passage over (or possibly tunneling through) the surface barrier.^{14} An important distinction between these phenomenological models and the more fundamental treatments of Kane,^{15} Mahan,^{16,17} and others is the assumption of an isotropic distribution of photoexcited electrons and isotropic scattering and the assumption that inelastic scattering depends only on the excited electron's energy.^{12} Such assumptions are not strictly correct, but as a practical matter, threestep^{13,18–23} and Momentsbased models^{9,11,24,25} have proven themselves expeditious and relatively accurate methods of modeling photoemission which are readily adapted to the needs of beam optics codes.^{4,26,27}
The moments $ M n ( k j ) = \u27e8 k j n \u27e9 $ can be used to find quantum efficiency QE and emittance $ \epsilon n , rms $ by utilizing the expectation values of the forward momentum $ \u210f \u27e8 k z \u27e9 $ or the mean transverse energy $ \u210f 2 \u27e8 k \u22a5 2 \u27e9 / 2 m $: for example, the quantum efficiency is
where $ R ( \omega ) $ is the reflectivity, D(E) is the transmission probability, and $ f \lambda ( cos \u2009 \theta , p ( E ) ) $ is the scattering factor governing losses given by
with $ p = \delta / l , \u2009 \delta ( \omega ) $ being the laser penetration depth, and $ l ( E ) = v ( E ) \tau ( E ) $ is the mean free path for an electron with a velocity v(E). In a similar vein, the emittance is
where $ \u210f k \u22a5 $ is the photonaugmented momentum transverse to the surface normal in vacuum^{28} and ρ_{c} is the radius of the illumination area. Other factors have their usual meaning.
In what follows, the units and terms are dictated by how the image charge barrier V(x) at the surface is represented, specifically
which follows the conventions established in Refs. 4, 29, and 30. Observe that the product of electric field $E$ and electron charge magnitude q is joined in the force term $ F = q  E  $, and similarly, q is combined with potential $\phi $ to give potential energy $ V ( x ) = q \phi ( x ) $. The image charge term $ \u2212 Q / x $ identifies $ Q \u2261 q 2 / 16 \pi \epsilon 0 = 0.36 $ eV nm. μ and $\Phi $ are the Fermi energy and work function, respectively. In the presence of a field, the maximum height of V(x) is reduced by the Schottky factor $ 4 Q F $ so that the “effective” work function is $ \varphi \u2261 \Phi \u2212 4 Q F $.
For metals, the Fermi Dirac distribution functions for low temperature become restrictions on the magnitude k such that the initial energy $ E = \u210f 2 k 2 / 2 m $ is restricted by $ E \u2264 \mu $ and $ E \u2265 \mu \u2212 \u210f \omega $. The transmission probability sets a restriction on the magnitude of the “forward energy” $ E x = E \u2032 \u2009 cos 2 \theta $ (with $ E \u2032 = E + \u210f \omega $ being the photoexcited energy) to be greater than the barrier height $ \mu + \varphi $. Semiconductors follow a similar development, but the electron affinity E_{a} and bandgap E_{g} usurp the roles of μ and $\varphi $ and induce changes in how E is defined. Under such approximations, for example, QE for metals is
where θ_{m} is the angle at which $ ( E + \u210f \omega ) \u2009 cos 2 \u2009 \theta m = \mu + \varphi , \u210f 2 k F 2 / 2 m = \mu $, and $ \u210f 2 k m 2 / 2 m $ defines k for which $ E + \u210f \omega = \mu + \varphi $. For semiconductors, Eq. (5) becomes
where $ \eta = cos \u2009 \theta $, with $ \eta m = cos \u2009 \theta m $ so as to show an alternate common formulation.
The isotropic and inelastic approximations allow treating delayed emission effects when the pulse duration of the drive laser becomes so brief that the depth of photoexcitation, the contribution of scattered electrons, and the magnitude of the relaxation time of the scattering processes^{10,31} matter in shaping the emitted bunch. For a flat surface, translational invariance means all photoexcitation and emission are independent of the y and z coordinates and dependent only upon the depth x of photoexcitation. Therefore, all photoexcitation at x can be modeled as coming from a point. Delayed emission is then modeled as shells of charge expanding radially outward from an origin until such time as they undergo a scattering event or pass the surface plane, as in Fig. 1. The theory and numerical implementation of that model is given in Sec. II.
II. MODELING DELAY FROM SCATTERED ELECTRONS
A delay model makes use of the following approximations and terms:

the initial energy of a photoexcited electron is (for metals) $ E = \mu + \u210f \omega $ or (for semiconductors) $ \u210f \omega \u2212 E g $, and therefore, the magnitude of the initial velocity v_{o} is the same for all;

time advances in increments with $ t n = n \Delta t $ denoting the present time and a time $ j \Delta t $ units earlier denoted by $ t n \u2212 j = ( n \u2212 j ) \Delta t $;

photoexcitations occur at fixed sites denoted by $ x l = \u2212 ( j \u2212 1 / 2 ) \Delta x $ for $ l \u2208 ( 1 \cdots N x ) $ with the first site a distance $ \Delta x / 2 $ from the surface; and

the number of times an electron scatters is s: a scattering event occurs every scattering time τ and the energy of the excited electron reduces by $ \Delta E $, resulting in a new velocity v_{s}.
Observe that Approx. (1) is tantamount to assuming that the barrier height is close to $ \varphi \u2248 \u210f \omega $ (metals) or $ E a + E g \u2248 \u210f \omega $ (semiconductors): the approximation is kept even as the difference between photon energy and barrier height increases. The meaning of some of the terms is illustrated in Figs. 1 and 2. The contribution of the l^{th} hollow shell (all charge is at the same radius from x_{l}) centered at x_{l} to the photoemitted current at a time t_{n} due to a shell created at a time t_{j} before it is found by summing up all of the current elements from all the created expanding hollow shells that satisfy emission conditions and that have undergone s scattering events. The elements of the calculation are codified in
where each term accounts for a specific process or feature, with the parameters described by the following:

J_{o} contains all of the various dimensioned factors grouped into a coefficient;

$ W n \u2212 j = I \lambda ( t n \u2212 j ) / I o $ is the ratio of the laser intensity at the time $ t n \u2212 j $ with a reference intensity I_{o};

$ S l s $ is a weighting factor specifying the shell charge in the l^{th} shell after s scattering events;

$ F l j s $ is the portion of the shell charge, or fractional charge, that would pass into vacuum if unobstructed due to the l^{th} shell at a time t_{j} from its creation and after s scattering events; and

$ D l j s ( n ) $ is the transmission probability at time t_{n} for electrons in the l^{th} shell created at a time t_{j} before the present time t_{n} for which the energy has been reduced by s scattering events.
The photoemitted current J_{total} at time t_{n} is then the sum over the matrix elements specified by $ J l j s ( n ) $. Each item in the list corresponds to a section below.
A. Coefficient
Let I_{o} be the maximum laser intensity of a pulse in the absence of noise. The number of absorbed photons in a time $ \Delta t $ per unit area for that intensity will be $ I o \Delta t / \u210f \omega $. The number of excited electrons in a time $ \Delta t $ per unit area will be $ ( J e \Delta t ) / q = ( I o \Delta t ) / \u210f \omega $. The amount of charge in the l^{th} shell will arise from the absorbed laser power at a depth x_{l} into the bulk, and its associated term will introduce a factor $ \Delta x $ into the numerator of the J_{o} coefficient that collects such factors. The shells expand with a velocity v_{o}, so the length scale for the fractional charge components will be $ v o \Delta t $. The spatial and temporal finite difference factors $ \Delta x $ and $ \Delta t $ will be related by the coefficient
The factor $ \Delta x / v o \Delta t $ can be rendered in terms of a ratio between the initial energy E_{o} of the photoexcited electron and a characteristic energy $ m \Delta x 2 / 2 \Delta t 2 $ by
B. Laser profile
The laser profile is $ W n \u2212 j \u2261 I \lambda [ ( n \u2212 j ) \Delta t ] / I o $. For a tophat pulse, $ W j = \Theta ( j ) \Theta ( N t \u2212 j ) $ for a pulse that began at time t = 0 and was of duration $ t pulse = N t \Delta t $, where $ \Theta ( x ) $ is the Heaviside step function. The Gaussian pulse is $ W j = exp ( \u2212 a ( j / N t ) 2 $, where a is set by a desired width of the pulse. An intermediate shape between the tophat and Gaussian pulse is of the form
where the coefficient a and the power p are chosen so that $ W N t $ is small and the pulse shape ranges from a tophatlike ( $ p \u226b 1 $) to a Gaussianlike ( $ p \u2248 2 $) profile. The behavior of Eq. (10) for various values of increasing p is shown in Fig. 3. The inclusion of noise, or “jitter,” is accomplished by creating a random number u for every W_{j} and modifying W_{j} according to
where $ 0 \u2264 u \u2264 1 $ and r_{n} is the maximum jitter.
C. Shell and fractional charge
The total photoexcited charge per unit time is $ \Delta Q = ( q / \u210f \omega ) I o \u2009 \Delta t $. Light absorbed as a function of depth into the bulk material, and so define
where x = 0 is the surface and $ x \u2192 \u2212 \u221e $ is deep into the bulk. The attenuation of light in the bulk material exponentially decays^{32} so that $ S 0 ( x ) $ is defined as
where $ \delta ( \omega ) $ is the laser penetration depth (e.g., 12.8 nm for copper) that governs the penetration of light into an attenuating medium, and the s = 0 superscript reflects that no scatterings have taken place. $ S 0 ( x ) $ is normalized so that its integral from $ \u2212 \u221e $ to 0 is 1. Discretizing by $ x \u2192 x l = \u2212 ( l \u2212 1 / 2 ) \Delta x $ therefore entails that the number of photoexcited electrons generated between $ x l \u2212 \Delta x / 2 = \u2212 l \Delta x $ and $ x l + \Delta x / 2 = \u2212 ( l \u2212 1 ) \Delta x $ is equal to
with $ r \u2261 \Delta x / \delta $ and so $ \u2211 l = 1 \u221e S l 0 = 1 $.
The 0superscript on $ S l 0 $ is an indication that no scatterings have occurred (s = 0) for the electrons making up the $ \Delta Q l $ of charge. The shell therefore expands with a velocity v_{s} and achieves a radius $ v s \tau $ when the s^{th} scattering occurs. To evaluate the next scattering event $ S l s $ given the previous event $ S l s \u2212 1 $, the contributions of all the expanding shells to the electron density at the various x_{l} must be summed, as indicated in Fig. 4: because of translational invariance along the $ y \u0302 $ and $ z \u0302 $ directions, scattering can be thought of as “resetting” the expanding shells, so that the shells collapse to points with a new weighting factor $ S l s $ and velocity v_{s} and begin expanding afresh, where the new $ S l s $ components are sums over the lsegments of the $ S s \u2212 1 $ shells (light blue and red vertical strips in Fig. 4 for two of the contributing shells). As the shell has a finite radius, define $ N s \u2261 v s \tau / \Delta x $, and so
where $ N \u2032 s $ excludes any points in vacuum (that is, $ l + k > 0 $ for all allowable k). The charge on the l^{th} segment is
after using $ r s \u2009 cos \u2009 \theta l \xb1 1 / 2 = x l \xb1 1 / 2 = \u2212 ( l \u2212 1 / 2 \u2213 1 / 2 ) \Delta x $, and so $ \Delta Q l $ is independent of the index l. Defining $ 2 N s + 1 = \Delta x / 2 r s $, Eq. (16) follows. Lastly, l + k must be positive because no shells are in vacuum, so that $ N \u2032 s < N s $ when l is small, that is, near the surface. The absence of those contributions causes $ S l s $ to appear cut off near the surface (most prominently for s = 1) and reflects the loss of electrons due to emission.
How much of the shell escapes into vacuum in a time $ \Delta t $ is the fractional charge $ F l j s $, where l marks the location of the center of the shell, j the time the shell has been expanding, and s the number of scattering events endured. If the shell passed into vacuum without impediment, then the fractional charge will be the charge on the circular ribbon with an arclength between the points $ v s t j $ and $ v s t j + 1 $ in Fig. 2. The current that passes x = 0 in a time $ \Delta t $ is $ \Delta Q l F l j s $.
$ F l j s $ is found from the difference between the amount of the shell in vacuum at time t_{j} and at time $ t j + 1 $. The fraction $ F ( \theta ) $ of a shell of radius r in vacuum as a function of angle θ is
Therefore,
Using $ cos \u2009 \theta j = x l / ( v s t j ) $, we get
where
and E_{s} is the energy of the electron after s scattering events.
For metals, electronelectron scattering dominates, whereas for semiconductors with a magic window, as in Fig. 5, electronphonon scattering dominates.^{13,31} For electrons originating at the Fermi energy μ (metal) or top of the valence band (semiconductor),
 Metal: The photoexcited electron can collide with an electron at the Fermi energy μ (max) or at $ \mu \u2212 E s $ (min). The average energy after collision is $ ( E s + \mu ) / 2 $ or μ, respectively, as in Fig. 5 (metal). Averaging the max and min cases gives$ E s + 1 = 1 2 { 1 2 ( E s + \mu ) + \mu } = \mu + \u210f \omega 4 s . $
 Semiconductor: phonon collisions reduce the energy to $ \Delta E $, as in Fig. 5 (semi), so that if E_{g} is the band gap,$ E s = E s \u2212 1 \u2212 \Delta E = \u210f \omega \u2212 E g \u2212 s \Delta E , $
where the final forms are found using mathematical induction. If $ f s \u2261 E s / E o $, then
where $ r \omega = \u210f \omega / ( \mu + \u210f \omega ) $ for metals and $ r \omega = \Delta E / ( \u210f \omega \u2212 E g ) $ for semiconductors. Therefore, $ v s = v o f s $.
D. Transmission probability
Photoemission fields are insufficient to cause quantum mechanical tunneling. The simplest approximation to the transmission probability is a step function of the form $ D ( E ) = \Theta ( E \u2212 V o ( t ) ) $, an instantaneous emission approximation that is valid so long as the transit time of the photoexcited electron over the barrier is much shorter than the fullwidth at half maximum (FWHM) duration of the laser pulse width. The transit time is compatible to $ v o \u2212 1 Q / F $ or less than a femtosecond and therefore much shorter than picosecondscale FWHM pulses. The approximation is therefore adequate. The time dependence of the barrier height is a consequence of the applied field F(t). For metals, using Eq. (4) gives the barrier height as
For semiconductors, the relation is more complicated because of the dielectric properties of the bulk material and changes to the nature of the emission barrier. The dielectric effects can be included by incorporating the dielectric constant K_{s} into the image charge via $ Q \u2192 ( K s \u2212 1 ) Q / ( K s + 1 ) $ as is commonly done.^{33} A more triangular barrier, however, induces changes in the transmission probability which depart from a stepfunction behavior.^{34} The energy in D(E) is $ E s \u2009 cos 2 \u2009 \theta $, where $ cos \u2009 \theta $ is the ratio of the distance to the surface x_{l} and the radius of the expanding shell $ v s t j $. That is, the argument of the Θfunction becomes $ ( 1 / 2 ) m ( x l / t j \u2032 ) 2 \u2212 \mu \u2212 \varphi ( t ) $, where $ j \u2032 \Delta t $ is the amount of time that the shell has been expanding since the last scattering event, and $ \varphi ( t ) = \Phi \u2212 ( 4 Q F ( t ) ) 1 / 2 $ for metals, with an analogous equation holding for semiconductors. Using the fact that if a and b are positive numbers, then $ \Theta ( a 2 \u2212 b 2 ) = \Theta ( a \u2212 b ) $, and then, the transmission probability element becomes
where $ t j = j \Delta t $ keeps track of how far back in time before the present time t_{n} is under consideration, and the denominator of the factor containing it arises from $ j \u2032 $.
E. Partial and total current
The current at time t_{n}, with no laser emission prior to t_{0}, due to all electrons that have scattered s times (the partial current) is found by summing over the spatial index l accounting for the contribution of all the shells in the bulk and the temporal index j accounting for the arrival of all shells beginning their expansion at a time t_{j} earlier than the present time t_{n}, or
where, in practice, the l summation extends only to N_{x} sufficiently large over the duration of the calculation, and $ S l s $ for $ l > N x $ is negligible. The total current is the sum over all scattering (partial) currents and is
where N_{g} is the number of scattering events that occur before all electrons have lost sufficient energy to no longer be emissioneligible: for a negative electron affinity (NEA) surface, N_{g} could be infinite, but for positive electron affinity (PEA) surfaces, it is always finite. If N_{g} is sufficiently large, then the scattered electron contribution acquires a diffusive character, modeled via an ad hoc means previously using a “Sphere” model:^{7} the formulation of the photoemitted current as arising from expanding shells of photoexcited charge enables a major correction to the previous model with its assumption of thermalized electrons and unvarying surface barrier. The consequences of a timevarying transmission probability on all populations that have scattered s > 0 times are now possible for the first time. The behavior of Fig. 6 creates an expectation borne out in simulation: the contribution of the partial term $ J s + 1 $ will in general be less than J^{s} because of the depletion of emissioneligible electrons near the vacuum interface. This occurs apart from energy losses due to scattering and a PEA barrier. The sum of all the s > 0 contributions when N_{g} is large results in behavior that is similar to the diffusive Sphere model, as shall be seen below, but it must be emphasized that the physics is decisively different.
III. DELAY MODEL RESULTS
A. Emission behavior
The partial currents $ J s ( n ) $ and total current J(n) are now evaluated. Three sets of evaluations are considered for metallike and semiconductorlike parameters as given in Table I. They are (i) the evaluation of the $ S l s $ factors for $ 1 \u2264 s \u2264 N g $ (Fig. 7); (ii) the evaluation of the $ J s ( n ) $ quantities for $ 1 \u2264 s \u2264 N g $ (Fig. 8); and (iii) the evaluation of the J(n) delayed emission quantities for timeindependent (F_{p} = 0) and timedependent ( $ F p > 0 $), as per Eq. (29) (Fig. 9). It is emphasized that the primary concern of the present study is how expanding shells distributed into the bulk material characterized by different expansion rates after any number of collisions nevertheless overlap and mutually contribute to the total current as governed by the relative magnitude of each scattering contribution contained in $ J s ( n ) $ [see Eqs. (27) and (28)] given that the shells encounter a barrier that dynamically changes. The effective mass variation, affecting as it does the magnitude of v_{s} and the scattering factor p in $ f \lambda ( \eta , p ) $, changes the weights of the terms in $ J s ( n ) $, whereas the primary concern in the theoretical development is instead how $ J s ( n ) $ is constructed from Eq. (7) in light of the timedependent processes, for which the magnitude of the effective mass emerges as a secondary concern. For purposes of illustration of the effects, therefore, the effective mass of the semiconductor was taken to be the same as the rest mass of the electron in vacuum so as to avoid wave function matching complications at the boundary, even though a smaller electron effective mass m_{n} entails a larger velocity and will lead to a longer delay time. This is not an unreasonable first approximation: effective mass contributions are more complex than simply assessing how $ ( k \u2212 1 \u2202 k E ( k ) ) $ behaves near the bottom of the conduction band. Preliminary results using density functional theory (DFT)^{35} indicate that in addition to a large scale free electron density of states, sharp features associated with extrema of energy bands arise and are characterized by a range of energies where electrons contributing to photoemission may be distributed as free electrons with a single particle effective mass. Such concerns, though, are outside the scope of the present study, for which establishing a lower bound is sufficient.
Parameter .  Symbol .  Unit .  Metal .  Semicon. . 

Wavelength  λ  nm  355  405 
Photon energy  $ \u210f \omega $  eV  3.493  3.061 
Eq. (10)  a  …  6  6 
Eq. (10)  p  …  32  32 
Eq. (29)  F_{o}  eV/μm  0  0 
Eq. (29)  F_{p}  eV/(μm fs)  100  20 
Jitter  r_{n}  …  10%  10% 
Dielectric const.  K_{s}  …  $\u221e$  6 
Work function  $\Phi $  eV  1.6  … 
Fermi energy  μ  eV  7  … 
Electron affinity  E_{a}  eV  …  0.4 
Band gap  E_{g}  eV  …  1.4 
Eq. (24)  r_{w}    0.3329  0.0602 
Relaxation time  τ  fs  2  16 
Laser penetration  δ  nm  12.6  40.0 
Spatial unit  $ \Delta x $  nm  0.25  1.4 
Temporal unit  $ \Delta t $  fs  0.25  1.4 
Spatial number  N_{x}  …  512  400 
Temporal number  N_{t}  …  512  400 
Eq. (9)  R_{o}  …  0.2709  1.3586 
Parameter .  Symbol .  Unit .  Metal .  Semicon. . 

Wavelength  λ  nm  355  405 
Photon energy  $ \u210f \omega $  eV  3.493  3.061 
Eq. (10)  a  …  6  6 
Eq. (10)  p  …  32  32 
Eq. (29)  F_{o}  eV/μm  0  0 
Eq. (29)  F_{p}  eV/(μm fs)  100  20 
Jitter  r_{n}  …  10%  10% 
Dielectric const.  K_{s}  …  $\u221e$  6 
Work function  $\Phi $  eV  1.6  … 
Fermi energy  μ  eV  7  … 
Electron affinity  E_{a}  eV  …  0.4 
Band gap  E_{g}  eV  …  1.4 
Eq. (24)  r_{w}    0.3329  0.0602 
Relaxation time  τ  fs  2  16 
Laser penetration  δ  nm  12.6  40.0 
Spatial unit  $ \Delta x $  nm  0.25  1.4 
Temporal unit  $ \Delta t $  fs  0.25  1.4 
Spatial number  N_{x}  …  512  400 
Temporal number  N_{t}  …  512  400 
Eq. (9)  R_{o}  …  0.2709  1.3586 
For semiconductors, the longer relaxation time of phonon scattering^{10,31} enables more $ S l s $ to contribute (larger values of s) and extend deeper into the bulk material. As more $ S l s $ contribute, the pulse shape of J(n) increasingly acquires a “driftdiffusionlike” characteristic reminiscent of similar tails observed in other simulations such as Monte Carlo^{36–38} or the Sphere model^{7} or experimentally (e.g., Fig. 3 of Ref. 39 and Fig. 9 of Ref. 5). Such tails have impact on the response time of photocathodes and affect their utility.^{40} Driftdiffusion components are particularly important for photocathodes with an NEA surface, such as GaAs,^{41} and contribute to their high quantum efficiency;^{13,37} analogous effects are seen in secondary emission in the reflection mode from hydrogenated diamond.^{42}
The contribution of the scattered components is most readily seen on comparing the profile of the incident laser pulse (with jitter) to the contributions of the partial currents $ J s ( n ) $ for each group of scattered electrons, as in Fig. 8. In the case of metals, only the very low work function (similar to a cesiated tungsten surface^{24}) allows any contribution to be seen for electrons that have scattered, although the approximation that all photoexcited electrons originated at the Fermi level likely overrepresents the population associated with the firstscattered electrons for metals. For semiconductors, where the energy loss per scattering is relatively small, many more events are required to remove excess energy (phonon collisions can also increase an electron's energy, but energy loss is the dominant mechanism, and the associated relaxation time reflects the cumulative behavior).
Lastly, delayed electrons experience a smaller barrier due to Schottky barrier lowering, and so, a time dependent simulation will reveal their contribution to the yield at later times. For simplicity, take the timevarying electric field to be linearly varying, or
The results of the simulation are shown in Fig. 9. For metals, the contribution is not significant, but for semiconductors, where a sizable fraction of the emitted electrons have undergone scattering, the effect is more noticeable. As a result, for applications requiring very fast pulses, the delayed emission effects can give a more sizable contribution: without a delayed emission model, understanding their impact on bunch shape in a beam optics code is not possible, even though the contribution of electrons emitted at inopportune times is a recognized issue.
The last observation becomes more complicated when models are modified by the presence of coatings, crystal structure, and changes to the emission barrier.^{8,35,43,44} Under such circumstances, resonances in emission can occur, analogous to (although different than) similar effects observed in field emission^{45,46} and discharges^{47} where the structure in the emission barrier dynamically modifies emission current. The presence of intentional heterostructures at the surface can therefore favor or suppress scattered electron contributions or affect the shape of the emitted pulse. For semiconductors, effective mass variations, if consequential, will change the emission probability to account for effective mass discontinuities when passing from material to vacuum and will therefore change QE and impact $ \epsilon n , r m s $. Such discontinuities will require modification of the matching equations in a Transfer Matrix Approach^{34,48,49} to evaluate D(E). When coupled with predictive modeling, these techniques can guide optimization of design parameters (including the film thickness and substrate composition). Given the recent focus on timeresolved measurements, strategies for achieving ultrafast pulses are a priority. The most straightforward method for increasing emission promptness is to reduce the cathode film thickness to be less than the optical penetration depth. For first order, this reduces QE (and, therefore, bunch charge), but techniques are emerging for optimizing emission from thin cathode films.^{50,51}
B. Performance in the beam optics code
Algorithms to evaluate the matrix elements of Eq. (7) and provide the timedependent emission were developed for the beam optics code MICHELLE. The electrostatic ParticleinCell (PIC) code MICHELLE^{26} has the following features: Image forces of electrons in the cavity walls are included. The selfB field of the beam is evaluated from Lorentz transformations. The rf fields within the cavities were determined using the AnalystMP eigenvalue solver^{52} and imported into MICHELLE. For the present simulations of a single cavity injector with a 2 mm radius spot size, the rf fields are imported into MICHELLE and a beam was launched based on a presumed laser pulse characterized by a tophat distribution in time and with uniform laser intensity across the cross section.
Two moments in the evolution of the bunch are contrasted: first, when the bunch is just emitted and second, just before the bunch enters the beam tunnel, as in Fig. 10 for the launching of a 1 nC bunch into an rf photoinjector modeled after a LANL normal conducting rf (NCRF) gun. The pulse shape is the same as those in the Shell/Sphere simulation of Ref. 7 to show the effects emission delay has. The beam is roughly barrel shaped at emission. The parameters are shown in Table II. A very similar behavior occurs in a different injector, where a sequence of images in Fig. 11 (the number in the lower right hand corner is the time index) again shows space charge distortion and tail formation.
Gun parameter .  Value .  Unit . 

Frequency  <1  GHz 
FWHM  22.4  ps 
Injection phase  $ 16.2 \xb0 $  … 
Peak $E$  $ \u223c 21 $  MV/m 
$E$ at $ 71 \xb0 $  $ \u223c 20 $  MV/m 
Laser energy  1.8  μJ 
Bunch charge  1.5  nC 
Gun parameter .  Value .  Unit . 

Frequency  <1  GHz 
FWHM  22.4  ps 
Injection phase  $ 16.2 \xb0 $  … 
Peak $E$  $ \u223c 21 $  MV/m 
$E$ at $ 71 \xb0 $  $ \u223c 20 $  MV/m 
Laser energy  1.8  μJ 
Bunch charge  1.5  nC 
Consider, then, the pulse shape and density profile of a bunch with and without the inclusion of delayed emission effects for two cases using semiconductor parameters (the effects are more pronounced) for parameters in Table III. Presented are the instances in a 2 D simulation just after emission (blue circle to left of Fig. 10) and just prior to entering the beam tunnel (red circle to right of Fig. 10). For the low bunch charge of 0.1 nC, the results just after emission are shown in Fig. 12, where the inclusion of the delay tail is clearly evident for the beam with the delay model employed. At the instant in time just as the beam is to enter the beam tunnel, as in Fig. 13, the beam shape and charge pattern due to the emission delay on the right show a prominent elongation and pulse distortion. The elongation is reminiscent of sphere contributions,^{7} but its abrupt truncation is new and indicative that only so many scattering events are tolerated before emission is suppressed. Moreover, the charge contained in the tail is no longer ad hoc but now predictively determined.
Gun parameter .  Value .  Unit . 

Frequency  <1  GHz 
FWHM  22.4  ps 
Injection phase  $ 16.2 \xb0 $  … 
Peak $E$  $ \u223c 10 $  MV/m 
$E$ at $ 16.2 \xb0 $  $ \u223c 2.9 $  MV/m 
Laser energy  0.15, 1.15  μJ 
Bunch charge  $ 0.1 , 1.0 0.1 , 1.0 $  nC 
Gun parameter .  Value .  Unit . 

Frequency  <1  GHz 
FWHM  22.4  ps 
Injection phase  $ 16.2 \xb0 $  … 
Peak $E$  $ \u223c 10 $  MV/m 
$E$ at $ 16.2 \xb0 $  $ \u223c 2.9 $  MV/m 
Laser energy  0.15, 1.15  μJ 
Bunch charge  $ 0.1 , 1.0 0.1 , 1.0 $  nC 
The behavior becomes exaggerated when the bunch charge increases by an order of magnitude to 1 nC, as in Fig. 14 for the beam just before it enters the beam tunnel: space charge forces clearly have a more dramatic impact beyond the greater content in the delayed part of the emission pulse. In these cases, the nominal laser pulse and resultant emitted current from the delay model are shown in Fig. 15.
IV. CONCLUSION
The present work develops and implements a substantial change to the theory of how delayed electrons contribute to a photoemission pulse by explicitly evaluating the current components of electrons that have scattered one or more times prior to emission. If the number of scatterings is sufficiently large, then the current from scattered electrons naturally acquires a diffusive character that can now be predictively evaluated rather than relying on an ad hoc diffusive current contribution (“Sphere” model) as in prior treatments. The new theory accounts for delays caused by, first, electron transport to the surface and, second, scattering changes to their kinetic energy: it thereby enables understanding the consequences of delayed emission effects utilizing beam optics codes.
When the laser penetration depth is sufficiently deep and the pulse width is sufficiently short (conditions that are of increasing technological interest), then even in the absence of scattering, delay effects can occur. With scattering, further delays are introduced, and the electrons undergo a change in energy depending on the scattering mechanism. Short pulse emission is of increasing technological interest for several reasons: high gradient injectors move to higher frequency which requires prompt psscale emission within the optimal rf injection phase; XFEL applications require aggressive bunch compression that motivates shorter pulses at the cathode; and emerging applications such as ultrafast electron diffraction and microscopy inherently require very short pulses (but much less bunch charge). While these applications have different tradeoffs, particularly between promptness and bunch charge, they both require a delayed emission model to inform optimization methods.
The model has been incorporated in the beam optics code MICHELLE and the effects of delayed emission quantified by the characterization of the emission tail and charges to the density distribution in the bunch. As a result, a teardroplike tail develops on the bunch although the tapered part of the drop is now truncated in contrast to the earlier Sphere model.^{7} When the surface field is timevarying, the delayed emission contributes to differences in the amount of emitted charge contributing to the beam bunch. As a result of the repartitioning of the charge density in the beam bunch both spatially and temporally from the emission delay, optimum system parameters such as timing of the laser pulse with respect to the electric field (strength and variation) both at the emitter surface and during the acceleration, the laser pulse shape, laser energy imparted to that target (affecting total charge emitted), and compensating magnetic fields may change in a significant way. Such variations may also affect beam head and tail evolution, as well as beam halos that may form, and having a predictive capability is important. The simulations have been performed using a conventional photoemission barrier (dielectricmodified image charge model): when the complexity of emission is altered by changes to the surface barrier model, then more dramatic departures will occur when delayed emission results are compared to conventional instantaneous emission models.
ACKNOWLEDGMENTS
The authors gratefully acknowledge support by the U.S. Department of Energy (DOE), Office of Science under the SBIR/STTR Grant No. DESC0013246.
References
Reference 9 applied the Moments approach to the evaluation of transverse emittance: by neglecting to include the contribution of photon energy to the transverse momentum of the emitted electron, the emittance was underpredicted by a factor of $\mu /(\u210f\omega +\mu )\u22480.77$ for copper parameters and 266 nm light. The neglected term was identified and corrected in Refs. 22 and 28.