Laser-cooled gases offer an alternative to tip-based methods for generating high-brightness ion beams for focused ion beam applications. These sources produce ions by photoionization of ultracold neutral atoms, where the narrow velocity distribution associated with microkelvin-level temperatures results in a very low emittance, high-brightness ion beam. In a magneto-optical trap-based ion source, the brightness is ultimately limited by the transport of cold neutral atoms, which restricts the current that can be extracted from the ion-generating volume. We explore the dynamics of this transport in a ^{7}Li magneto-optical trap ion source by performing time-dependent measurements of the depletion and refilling of the ionization volume in a pulsed source. An analytic microscopic model for the transport is developed, and this model shows excellent agreement with the measured results.

## I. INTRODUCTION

Photoionization of laser-cooled atoms offers a novel pathway for constructing high-brightness ion sources.^{1–7} These sources are fundamentally enabled by the microkelvin temperatures achievable through laser cooling, where atomic ensembles are produced with a very small momentum spread. Upon photoionization and extraction, the low momentum spread in the transverse direction results in a very low emittance ion beam. This emittance, when combined with a reasonable amount of current, results in a brightness that rivals or surpasses conventional ion sources.

Cold atom ion sources benefit from several advantages over conventional tip-based ion sources, such as access to new ionic species, inherent isotopic purity, insensitivity to contamination, and a low energy spread, which reduces chromatic aberration in ion optical systems and enables source operation at low accelerating voltages. Applications for these new sources are broad, encompassing milling, imaging, spectroscopy, and implantation. For example, a Li magneto-optical trap ion source (MOTIS) has demonstrated high surface sensitivity in imaging,^{8} imaging of optical modes in nanophotonic resonators,^{9,10} and utility in studying Li-ion-battery-relevant materials at the nanoscale.^{11,12} Cold atom ion sources of heavier species such as Cs enable nanometer scale milling resolution, and improved secondary ion mass spectroscopy.^{13} Recent work has demonstrated active feedback control of an ion beam using a measurement of the momentum of the photoemitted electron to correct the trajectory of each generated ion.^{14,15} Such control raises the possibility of deterministic implantation of single ions, as well as active compensation of the random transverse energy of the beam to further increase the resolution of these instruments.

Several methods exist for laser cooling atoms into cold atomic beams or trapped gases to realize an ion source. The technical complexity and utility of each approach varies with atomic species, as does the path toward optimizing ion source brightness. In the case of Cs, an atomic beam source with multistage laser cooling has been demonstrated with peak brightness as high as 2.4 × 10^{7} A m^{−2} sr^{−1} eV^{−1}, a value 24 times greater than the typical Ga^{+} liquid metal ion source brightness.^{13} Li, on the other hand, presents significant complexities for multistage laser cooling due to its high Doppler temperature and repumping requirements, and it is more practical to create and optimize a source on a simpler platform using photoionization in a magneto-optical trap (MOT). Although MOTIS brightness is fundamentally limited by transport of cold atoms into the ionization volume, reasonably high brightness can nevertheless be achieved with Li. For example, a Li MOTIS has been optimized to have a peak brightness of 1.2 × 10^{5} A m^{−2} sr^{−1} eV^{−1} and current as high as 1 nA (at lower brightness) using enhanced MOT loading and two-photon ionization.^{16}

Because of the range of applications that a Li MOTIS can address, it is useful to investigate in some detail the transport processes that limit the brightness. A full understanding of this process will not only enhance intuition for future MOT-based system design but also allow detailed modeling for optimizing beam brightness. To this end, we study here the time dependence of the depletion and refilling of the ionization volume in the Li MOTIS described in Ref. 16. We experimentally measure the current extracted in pulsed operation and compare these measurements with a simple hard-sphere model and an analytic solution to the transport equation. Related transport dynamics have been explored in Ref. 5 with emphasis on optimizing quasicontinuous brightness, considering long-distance transport in a MOT where cooling is shut off during ionization and the excitation laser beam has a significant pushing effect. Here, we concentrate on transport localized to the ionization region in a MOT that is not turned off during ionization and find that ballistic motion dominates on the length scale relevant to ionization even in the presence of active cooling forces.

## II. EXPERIMENTAL METHODS

We realize a lithium MOTIS by photoionizing neutral ^{7}Li atoms at the center of a 3D MOT, as described in Ref. 16. The MOT is continuously loaded from a 2D MOT with pusher beam^{17} fed by a Li evaporator, producing a slow (≈30 m/s) flux of atoms. Loading of the MOT occurs at a rate of ≈1.1 × 10^{9} s^{−1} (equivalent to an ion current of 180 pA). The atoms in the MOT have an estimated temperature of ≈950 *μ*K and an estimated peak density of ≈4.2 × 10^{17} m^{−3}. The background pressure in the vacuum chamber is ≈10^{−5} Pa (≈10^{−7} Torr) when the evaporator is turned on, resulting in a MOT lifetime of 240 ms. The total MOT population without ionization is ≈2.7 × 10^{8} atoms. The MOT has a nearly Gaussian density profile with 1/*e*^{2} radius 480 *μ*m (830 *μ*m) along the strong (weak) direction of the magnetic field gradient. The magnetic field for the MOT is created by two stacks of Nd_{2}Fe_{14}B ring magnets mounted outside the chamber, which create a quadrupole field with gradient ≈0.5 T m^{−1} (0.25 T m^{−1}) in the axial (radial) direction. The axis of the ring magnets is oriented perpendicular to the ionization direction.

The laser cooling and trapping for the MOT is done with two laser frequencies obtained from a single Ti:sapphire laser utilizing an acousto-optical modulator (AOM). The two laser frequency components are superimposed using a fiber combiner. One component is tuned five linewidths (5Γ) below the $F=2$ to $F\u2032=3$ D2 transition at 671 nm and has an intensity of 7.5 mW/cm^{2}. The other component, used to repump atoms that have decayed to the $F=1$ hyperfine state, has a frequency tuned 2Γ below the $F=1$ to $F\u2032=2$ D2 transition and has an intensity of 5.8 mW/cm^{2}. The laser light is split into three counter-propagating pairs of beams, which are circularly polarized according to standard MOT practice.^{18}

Ionization takes place in a two-photon process using two perpendicular laser beams, each tightly focused to the same spot. Ground state 2S, $F=2$ atoms are excited to a 3P_{3/2} state using a UV excitation beam with a wavelength of 323.36 nm. A second laser beam at ≈799 nm produces near-threshold photoionization of atoms in the 3P state. The overlap of these two focused beams creates a roughly spheroidal ionization volume. The UV excitation beam is focused to ≈6.7 *μ*m (1/*e*^{2} radius) with ≈17.4 *μ*W of power, and the IR ionization beam is focused to a 1/*e*^{2} radius of ≈9.2 *μ*m with ≈380 mW of power. Ions formed at the center of the MOT are accelerated in an electric field to create an ion beam, as shown in Fig. 1. The MOT is located between two electrodes, one at +3 kV and the other held at ground to create a ≈44 kV/m electric field at the atoms (note the DC Stark shift due to this field is negligible and does not affect the behavior of the MOT). The accelerated ions travel through an einzel lens and are either imaged on a multichannel plate (MCP) located ≈38 cm downstream or collected in a movable Faraday cup.

The transport dynamics in the MOT are observed by using an AOM to modulate the UV excitation light with 1/*e* rise time <1 *μ*s, effectively switching on or off the ionization rate rapidly compared to changes occurring in the local MOT population. The time dependence of the ion current is observed by switching the UV light on for 30 *μ*s while simultaneously ramping a deflection voltage to create a streak camera image on the MCP. The deflection is achieved by moving the Faraday cup to the side and applying the ramped voltage to its outer electrode. This causes the ion beam to sweep across the MCP surface, creating a spatial map of the time dependence of the ion current. A complementary metal-oxide semiconductor (CMOS) camera records the MCP's phosphor screen, and the resulting image is corrected for spatial variation in detection gain across the MCP. The einzel lens is adjusted to focus the ion beam onto the MCP in order to minimize the spot size and the corresponding temporal resolution. With a deflector scan rate of 130 kV/s, the ion beam spot size results in a temporal resolution of 1.4 *μ*s (one standard deviation of a Gaussian fit). Scanning of the deflector is synchronized with the switching of the UV excitation light to enable integration over ≈500 repetitions for a single camera exposure. The repetition period is set at 1.5 ms to ensure sufficient time for the MOT to attain steady state in between ionization sequences. In addition to observing the ion beam's time dependence, we also measure the average current with the Faraday cup positioned to capture the entire beam, using an ammeter with a time constant much longer than any switching of the ionization rate. As described below, this provides a quantitative measure of the MOT refill dynamics when a suitable pulse sequence is applied to the ionization.

## III. EXPERIMENTAL RESULTS

In order to focus on cold atom transport limitations, we have chosen to operate our Li MOTIS in a regime where such transport effects are most prominent. To reach this regime, we use tightly focused ionization lasers with high powers to create a high ionization rate and use pulsed ionization over time periods that are short compared with the MOT lifetime to avoid complications due to changes in overall MOT population. In this regime, we hypothesize that the MOTIS behavior is as follows.

When the ionization volume at the center of the MOT is first exposed to a UV pulse in combination with the (continuously on) IR beam, the instantaneous extracted current is proportional to the ionization probability and the unperturbed local MOT density. Very quickly, though, ionization causes density depletion in the ionization volume and surrounding region because cold atom transport cannot replenish the density quickly enough to match the ionization rate. The current then decreases as the density depletes. After some time, a new quasiequilibrium is reached when the reduced ionization rate matches the cold atom transport rate. When the pulse turns off, the depleted ionization volume begins to refill due to cold atom transport, eventually reaching the original unperturbed state of the MOT.

Using the streak camera method described above, we have measured the time dependence of the current during a two-pulse sequence to confirm this hypothesis. Figure 2 shows the results of this measurement. As expected, the pulses show a large initial transient current followed by a lower, quasi-steady-state current. Since the magnitude of the transient response provides a measure of the density within the ionization volume at the beginning of each pulse, we use the height of the transient response in the second pulse to measure the recovery of the density depletion after the first pulse ends. Repeating the experiment for varying delay times between the first and second pulses, we obtain the time dependence of the density depletion refill by plotting the ratio of the second to the first peak transient response versus delay time. The resulting measurements are shown in Fig. 3 (red squares).

Also shown in Fig. 3 are the results of the time-averaged current measurements using the Faraday cup, conducted with the same two-pulse sequence (black circles). For each value of the delay between the two pulses, we collect and average the current over 3000 repetitions of the sequence, repeating at an interval of 1.5 ms. We then determine the average current during the 60 *μ*s that the ionization is on by correcting for the 4% duty cycle. As with the peak ratio data, the average current has a smaller value at short delay times and then increases as the delay time is increased until it reaches an asymptotic value. This behavior is consistent with our model of the transport dynamics: the first pulse depletes the ionization volume and then the second pulse is smaller if it arrives before the volume has a chance to refill, resulting in a smaller average current. If the delay time is long enough, then the volume refills to the unperturbed MOT density by the time the second pulse arrives, and the average current becomes independent of delay time.

## IV. THEORY

Given the clear measurements of ionization volume refill we have obtained, it is instructive to develop a simple theoretical description of the process. We begin with an ansatz where the MOT population consists of a cloud of cold, noninteracting atoms that is essentially uniform over the region of interest around the ionization volume. We further assume that the motion of the atoms is purely ballistic, the temperature remains constant and uniform, and all magnetic restoring forces are negligible because the MOT is located close to the magnetic field zero. In reality, atoms in a MOT are embedded in an optical molasses and so should be considered as moving in a diffusive medium. However, the diffusion constant is such that the mean free path is much larger than the size of the ionization volume, so ballistic transport is the more appropriate picture for our situation (see Appendix A). The ionization process is treated as a pure loss channel with an effective rate, and no attempt is made to model the details of the optical transitions either in the ionization or in the MOT laser cooling.

### A. Hard-sphere model

A certain amount of intuition into the time dependence of the refill dynamics in the MOT can be gained by considering a “hard-sphere” model, in which the ionization is assumed to occur entirely within a sphere of radius $a0$. In this model, we can write a simple differential equation for the number of atoms in the ionization volume, $N(t)$,

where $N0=4\pi 3a03n0$ is the initial population, $\beta =3v\xaf/(4a0)$ with $v\xaf=8kBT/\pi m$, and $\alpha $ is the ionization rate. Here, $n0$ is the unperturbed MOT density and $v\xaf$ is the most probable atomic velocity [$kB$ is Boltzmann's constant, *T* is the MOT temperature (950 *μ*K), and *m* is the atomic mass]. The three terms on the right-hand side of Eq. (1) are the ballistic flow into the sphere, the ballistic flow out of the sphere, and the loss due to ionization, respectively, where we have used the formula $14nv\xafA$ for the flow across a surface area *A* in an ideal gas of density *n*.^{19}

Equation (1) can be solved explicitly for the ionization pulse sequence used in the experiment with an initial pulse of length $\tau p$, a refill period of length $\Delta \tau $, and a second pulse of length $\tau p$. In so doing, we find that both the ratio of the transient peaks and the average Faraday cup current take the form $A0\u2212A1e\u2212\beta \Delta \tau $, where $A0$ and $A1$ are large, complicated functions of $\alpha ,\beta ,\tau p$, and $N0$. This solution isolates the dependence on $\Delta \tau $ into a simple exponential and provides a simple model with which to fit the data in Fig. 3. As seen in the figure, the fit is excellent and setting the time constant of the fit equal to $\beta \u22121$ results in *a*_{0} = (14.6 ± 0.4) *μ*m for the peak ratio data or *a*_{0} = (15.3 ± 0.5) *μ*m for the Faraday cup data. These two numbers, in agreement within the uncertainty, correspond to a reasonable radius that is commensurate with the beam waists of the ionization lasers (6.7 *μ*m and 9.2 *μ*m $1/e2$ radius for the UV and IR lasers, respectively).

### B. Boltzmann equation

While the hard-sphere model provides some intuition into the time dependence of the refill, its treatment of the transport oversimplifies the spatial dependence of the density and ionization probability. As a result, it does not do a good job modeling the coefficients *A*_{0} and $A1$. These coefficients are actually not independent in the model, as they both depend on the basic parameters $a0,\alpha ,n0,v\xaf$, and $\tau p.$ Among these parameters, only $a0$ and $\alpha $ can be considered as free fitting parameters for the given experimental conditions since $n0,v\xaf$, and $\tau p$ are measured quantities. With only two free parameters, a good fit to the data in Fig. 3 is not possible, as shown by the dashed line in the figure.

To more accurately model the spatial dependence of the density and the ionization probability, we use a collisionless Boltzmann equation with a loss term to solve for the distribution of atoms in the MOT as a function of time. For the ionization process, we assume a spherically symmetric Gaussian ionization profile. Assuming spherical symmetry is a somewhat crude approximation to the crossed laser beam intersection (see the actual observed current distribution shown in Fig. 1) but it nevertheless contains the essential spatial scale of the ionization and allows reducing the problem to a one-dimensional solution. We write the equation governing the distribution $f(r,v,t)$ of atoms in the MOT as

where $r$ is the spatial coordinate, $v$ is the velocity, *k* represents the peak ionization probability per unit time, and $\sigma $is the standard deviation of the ionization profile. Because this equation does not have any explicit derivatives in $v$, it can be solved for a specific $v$ and then averaged over all possible values of $v$, assuming the velocity distribution is an unperturbed Maxwell–Boltzmann distribution. Details of the solution are given in Appendix B. Assuming the ionization pulse turns on at $t=0,$ the spatial density distribution during the pulse can be written as

where $f0$ is a normalization constant and $vth$ is the thermal velocity $kBT/m$. Note we use spherical coordinates in both position and velocity space to take advantage of the (assumed) spherical symmetry, and the normalization factor $f0$ includes the concomitant factor of $4\pi .$ To obtain the total current as a function of time, we multiply the density by the ion charge *q* and the Gaussian ionization probability and integrate over space,

After the pulse is over and the depleted region begins to refill, the distribution evolves according to Eq. (2) with $k=0$. The solution in this case is simply $f(r\u2212vt,v,0)$, where $f(r,v,0)$ is the distribution at time $t=0.$ Taking the distribution at the end of the pulse as the starting distribution, we write

where $\tau p$ is the pulse duration. When the second pulse arrives, the distribution again evolves according to Eq. (2) with $k\u22600$ but with an initial distribution given by Eq. (5) with $t=\Delta \tau $. Now, the density as a function of time is

where *t* is measured from the start of the second pulse. Carrying out the integrations in Eqs. (3)–(6) numerically, we plot the results in Fig. 4, using the values for *k* and $\sigma $ determined from a fit to the electrometer refill data (see below) and $vth=1.06m/s$, corresponding to $T=950\mu K$. Figure 4(a) shows the time dependence of the current during the first 30 *μ*s pulse, where we have convolved the theory with a Gaussian function with standard deviation 1.7 *μ*s to account for combined streak camera and AOM resolution and normalized the peak of the transient to a value of 1. Comparing the theoretical time dependence to the streak camera data shown in Fig. 2, we see that the theory reproduces qualitatively the expected behavior: an initial peak in the current, followed by a drop off to about 70% of the peak value. A qualitative comparison with the experiment is all that can be done here because the time offset and resolution line shape cannot be accounted for satisfactorily. Figure 4(b) shows the relative density during the initial pulse as a function of radial coordinate *r* at $t=(0,1,2,5,10,and30)\mu s$. This plot shows how the depletion starts with a profile reflecting the ionization profile but grows deeper and wider as ionization continues. Figure 4(c) shows the relative density during the refill period at delay intervals $\Delta \tau $ of $(0,1,2,5,10,and30)\mu s$ after the first pulse has finished. Here, we see a gradual refilling of the distribution, as expected. Depletion of the distribution during the second pulse is very similar to the depletion during the first pulse, albeit starting from a lower initial level when $\Delta \tau $ is small, and so is not shown in Fig. (4).

To calculate the average current as would be seen by the electrometer, we write

and calculate the total charge in both pulses as

The average current is then $Iav(\Delta \tau )=Q(\Delta \tau )/2\tau p$. $Iav$ can be numerically evaluated with *k*, $\sigma $, and $f0$as parameters. Allowing these three parameters to be free, it is now possible to obtain a good fit of $Iav$ to the electrometer refill data, and the results are shown in Fig. 5. We note that the scale factor $f0$ could in principle be fixed though normalization but that would require good knowledge of the peak MOT density at the center, and also any detection inefficiencies in the ion current measurement, both of which we do not have good knowledge of. For this fit, *k *= (0.164 ± 0.032) *μ*s^{−1} and *σ* = (11.2 ± 0.8) *μ*m (uncertainty is one standard deviation for the fit). As can be seen, the theoretical curve matches the experiment extremely well. The optimized value of $\sigma $ is somewhat larger than the UV and IR beam radii (6.7 *μ*m and 9.2 *μ*m, respectively) but considering that the UV transition is highly saturated, and the actual beam overlap is more complex than a pure spherically symmetric Gaussian, the value is quite reasonable. The optimized value of the ionization rate$k$ is also quite reasonable considering the following. The total current can be estimated by noting that it should be equal to $k(2\pi )3/2\sigma 3n0,$ where $n0$ is the MOT density at the center. Combining the optimized values for *k* and $\sigma $ with $n0=3.3\xd71017m\u22123$ (a value lower than the nominal density because of depletion) results in a total current of 190 pA, in agreement with the magnitude of the current shown in Fig. 5.

## V. SUMMARY AND CONCLUSIONS

We have presented measurements of the transport dynamics within a Li MOTIS in the form of the time dependence of the current during pulsed ionization and the refilling of the depleted ionization volume between two pulses. We have modeled these dynamics via simple ballistic transport of cold atoms and obtained an analytic solution that agrees extremely well with the measurements when three basic parameters are adjusted to match theory with experiment. The optimized values of these quantities are within expectations for the experiment.

By demonstrating good agreement between a simple ballistic theory and observed behavior of the MOTIS, we have established that ballistic transport is a suitable approach toward analyzing the microscopic behavior of a Li MOTIS, understanding its limitations, and optimizing its performance. This result is important for future development and optimization of cold atom ion sources, as it establishes a pathway toward a more complete simulation framework. Having a good understating of the cold atom dynamics in these sources is crucial input for such a framework, and when it is integrated with existing ion optical models, a complete picture of source performance will emerge, enabling innovative designs.

While the present work gives support to a simple ballistic transport model, it is important to recognize that we have only seen this so far for the limited performance range of the current experimental parameters in a Li MOTIS. Magneto-optical traps have a wide dynamic range of operation, and their behavior differs significantly depending on atomic species, laser power and detuning, magnetic field gradient, lifetime, and geometry. This is particularly true if auxiliary cooling schemes are used, such as dark state cooling,^{20} gray molasses,^{21} or cavity cooling.^{22} The possibility of diffusion-limited dynamics in the presence of laser cooling cannot be ignored for all situations and may play a significant role in systems where extremely low temperatures are sought. In some cases, neither diffusion nor ballistic transport will dominate, and it will be necessary to solve a more complicated set of equations to model cold atom source behavior. In any case, our work here represents a step in the progression toward full modeling of cold atom source behavior and will hopefully stimulate further work in this area.

## ACKNOWLEDGMENTS

J.R.G. acknowledges support from the National Research Council's Research Associateship Program and W.R.M. acknowledges support under the Cooperative Research Agreement between the University of Maryland and the National Institute of Standards and Technology Physical Measurement Laboratory, Award No. 70NANB14H209, through the University of Maryland.

### APPENDIX A

We present here a justification for assuming ballistic transport rather than diffusion for describing the transport dynamics in our experimental situation. While the actual optical molasses diffusion constant within the MOT is difficult to estimate, we can argue as follows: if the current in a MOTIS were limited by diffusion with constant *D*, the steady-state density *n* near the center of the MOT would obey the equation $D\u22072n(r)\u2212L(r)=0$, where $L(r)$ is the spatially dependent loss due to ionization. This equation can be solved assuming delta-function ionization at the origin, $L(r)=I\delta (r),$ where *I* is the total current (in ions per second). The result is

where we have used the boundary condition that the density far from the ionization region is $n0$, the nominal unperturbed central density of the MOT. Because $n(r)$ cannot be negative, Eq. (A1) would imply that $D>I/(4\pi n0r)$ for any *r*. Inserting typical observed values $I=100pA$, $n0=4\xd71017m\u22123$, and $r=10\mu m$, we obtain $D>1.2\xd710\u22125m2/s$. However, this value of the diffusion coefficient results in a mean free path of $\lambda =2D/vth=24\mu m$, where $vth=1.06m/s$ is the thermal velocity at a temperature of 950 *μ*K. Since diffusion-limited motion assumes the mean free path is much shorter than the spatial dimensions of interest, the fact that $\lambda >r$ indicates diffusion is not a good description of the process, and collision-free ballistic motion is a more accurate assumption. A further confirmation of this long mean free path can be obtained from Li MOT measurements in the literature. Schünemann *et al.*^{23} measured a relaxation time of 100 *μ*s in a MOT with radius 130 *μ*m and temperature 300 *μ*K. These can be combined to estimate a diffusion constant^{24} of $1.7\xd710\u22124m2/s$ with a corresponding mean free path of 570 *μ*m.

### APPENDIX B

We describe here the method used to solve Eq. (2), the collisionless Boltzmann equation with a spherically symmetric Gaussian loss term. Because Eq. (2) has no explicit derivatives in $v$, we can treat $v$ as a parameter and assume the velocity distribution remains unchanged as a Maxwell–Boltzmann distribution with standard deviation $vth=kBT/m$. We write Eq. (2) in cylindrical coordinates and take the $z^$ axis to be along **v**,

This is a one-dimensional time-dependent differential equation in *z*, which has a general solution,

with $G(\rho ,z\u2212vt,v)$ being any arbitrary function of $\rho ,$$z\u2212vt$, and *v* chosen to yield the desired boundary condition at $t=0$. For a solution during the first laser pulse, we choose $f(\rho ,z,v,0)=f0e\u2212v2/2vth2$, i.e., a uniform distribution in space and Maxwell–Boltzmann distribution in velocity ($f0$ is a normalization factor), which results in

In order to determine the quantities of interest, such as the spatial distribution or the total current as function of time, we need to average over the velocity. To do this, we convert Eq. (12) to spherical coordinates via $\rho =rsin\theta $ and $z=rcos\theta $ and write

For a solution after the first laser pulse is finished at time $\tau p$, we desire a solution to Eq. (B1) with $k=0$. This is simply any function of $z\u2212vt$, where *t* is measured from the end of the pulse. We use $frefill(\rho ,z,v,\tau p,t)=f1(\rho ,z\u2212vt,v,\tau p)$ so that the distribution is correct at $t=0.$ For the distribution at time *t* after the second pulse arrives, we use Eq. (B2) again, this time requiring that the distribution at $t=0$matches the distribution at the end of the delay, i.e., $f2(\rho ,z,v,\tau p,\Delta \tau ,0)=frefill(\rho ,z,v,\tau p,\Delta \tau )$. We write

The distributions $f1(\rho ,z,v,t),frefill(\rho ,z,v,\tau p,t)$, and $f2(\rho ,z,v,\tau p,\Delta \tau ,t)$, corresponding to the periods during the first pulse, between pulses, and during the second pulse, respectively, can be converted to spherical coordinates and integrated over velocities as was done for Eq. (B4) to generate the spatial distributions as a function of time in these periods. These spatial distributions can also be integrated over space to provide the total current as a function of time.