2D Planar PIC Simulation of Space Charge Limited Current With Geometrical Parameters, Varying Temporal-Profile and Initial Velocities

In this paper, we simulated the space-charge effect current by the full electromagnetic particle-in-cell (PIC) code UNIPIC in a planar gap. Compared with Lau’s numerical work, we introduce the applied voltage <inline-formula> <tex-math notation="LaTeX">$U$ </tex-math></inline-formula>, gap separation <inline-formula> <tex-math notation="LaTeX">$D$ </tex-math></inline-formula>, suppression magnetic field <inline-formula> <tex-math notation="LaTeX">$B$ </tex-math></inline-formula>, electron-injection energy spectra <inline-formula> <tex-math notation="LaTeX">$f(E)$ </tex-math></inline-formula> and temporal profiles <inline-formula> <tex-math notation="LaTeX">$f(t)$ </tex-math></inline-formula> as the factors affecting the current densities and the particle distributions. It is shown that <inline-formula> <tex-math notation="LaTeX">$U$ </tex-math></inline-formula> and <inline-formula> <tex-math notation="LaTeX">$D$ </tex-math></inline-formula> could change the particle distribution significantly, causing the electrons’ escape, while they bring little impact on the ratio of 2-D Child-Langmuir current to the classical 1-D value <inline-formula> <tex-math notation="LaTeX">$J_{\mathrm {CL}}(2)/J_{\mathrm {CL}}$ </tex-math></inline-formula> <xref rid="deqn1" ref-type="disp-formula">(1)</xref>. Considering the temporal profile, it is interesting that the transported current density is inversely going to decrease when the injection rise time is much smaller than the transit time, which explores the smaller temporal profile numerically compared with Valfells’s experiments. Last, the preliminary simulation of the combined effect of temporal profile and initial electron energy spectra is presented, and it is evident that the thermal electrons obeying the Boltzmann distribution generate a larger transported current density than the photoelectron under exponential distribution, which might be an inspiration to the engineering application.


I. INTRODUCTION
Owing to the force of intense electric field, thermal effect or photons knocking [1]- [4], electrons will emit from the cathode, forming the electron flow under the externally applied voltage. Although this problem has been studied for almost over a century, it is still a subject of much interest in electron guns, accelerator physics, sheath physics, vacuum microelectronics, and high-power diodes [5]- [7]. Because of the complex interaction between charged particles and the external fields or self-consistent field, the emitted electrons usually maintain a dynamic state [8]- [10]. As the number of emission electrons grows large, the resulting electric field is significant to reflect the injected electrons, limiting the maximal transported current which is also called SCL current. Then the virtual cathode [11] generates when the quantity of particles continues increasing, giving rise to the virtual The associate editor coordinating the review of this manuscript and approving it for publication was Giovanni Angiulli . cathode oscillation which is widely applied to generate high power microwaves [12]- [14].
The famous three-halves power Child-Langmuir law gives the one-dimensional (1-D) time-independent SCL current density by two integrations of Poisson's equation in an infinite planar disk with gap separation D and gap voltage U [15], [16] where e, m, and ε 0 are the electron charge, electron mass, and the free space permittivity, respectively. In the case of classical 1-D SCL emission, the normal electric field at the cathode is set zero due to plenty of electrons reserving back to the emitter surface. However, the zero-potential gradient condition is invalid under the circumstance of electrons with finite initial velocity. Hence, the assumption that the potential gradient at the cathode is zero can only solve the problem involving electrons emitted with zero initial velocity [1].
While the actual initial velocity should be related to the particular physical emission mechanism. Jaffé investigated VOLUME 10, 2022 This work is licensed under a Creative Commons Attribution 4.0 License. For more information, see https://creativecommons.org/licenses/by/4.0/ the 1-D effect of electrons with uniform initial velocity on the SCL current density several decades ago [17], and Liu et al. recently studied the initial velocity effect on SCL current under the approximation that the initial injection velocity is independent of the emission mechanism and is represented by an average velocity v 0 . Nevertheless, the distributions of electron injection velocities affect the field and current in the diode a lot [9]. In this case, our team has summarized a new SCL current formula to display the varieties of J CL (2)/J CL (1) with initial velocities described by Maxwell distribution [18].
More and more researchers focus on the 2-D geometrical effect in the last two decades, and not only analytic descriptions [19], [20] but also numerical simulations [21]- [24] and experiments [25] have been reported. In these researches, the 2-D geometrical effect has been proved to increase the SCL current density, and many different empirical formulas have been concluded to describe the relation of J CL (2)/J CL (1) for strip model and disk model [21]- [24]. Additionally, there have been other analytical and simulation works on figuring out the effects of initial velocity, limited emission area, finite but varying pulse duration, external magnetic field, bipolar flow, and even the quantum effects in 2-D diode gap [1], [11], [21], [25]- [30]. For the temporal effects, Griswold et al. provides the upper bound of the SCL current [31] and Zhu et al. gives the current density limit in the Coulomb blockage regime [32]. Besides the current in vacuum, a simplified theory of relativistic space-chargelimited current in 2D thin-film massive Dirac fermions has been formulated [33]. Recently, the three-dimensional studies have been carried out to discover the effect of the microfabricated structures near the emitter in the classical and asymmetrical planar diode [34], [35] as well as some overall effects of nano emitter on the SCL current [36], [37].
Although most of the diode parameters have been studied to obtain their effects on the SCL, there still are some issues unsolved to be investigated. As is known to all, Lau et al. have established a comparably complete system in this area and given the law of the geometrical factor [19], [21], such as W /D for strip model and R/D for disk model, which will bring great change to the SCL current density in the 2-D model due to the edge-effect. Nevertheless, no data have demonstrated the effects of other parameters, such as the gap separation, applied electric or magnetic field can be ignored and the universality of the scaling law has not been cleared. Hence, this paper simulates the particle distribution and SCL current density with different parameters under different fixed geometrical factor W /D for a strip diode by the particle-in-cell (PIC) method. The combinational effect of temporal profile and emission velocity on SCL current in the strip model are computed as well. Section II briefly introduces the computational model. In Section III, the effects of voltage and gap separation, and the particle distributions are given. Section IV presents the combinational effect of temporal profile and different electron injection spectra on SCL current. Section V summarizes this work.

II. SIMULATION MODEL
PIC method is used here to study the SCL current. Poisson's equation Eq. (2) and Newton-Lorentz equation Eq. (3), solving the electromagnetic fields and the motion of charged particles, are the basic equations, where , E, ρ and ε are the potential, electric field, electric charge density, and permittivity, respectively. Since our work needs to consider the effect of initial velocity and applied voltage, the energy of particles might exceed 2.57 keV, which means their velocity could be over one tenth of light speed. So, Newton-Lorentz equation Eq. (3) should be in relativistic form, where r, v, m, q, and B are displacement, velocity, electron mass, electron charge and magnetic field, respectively, and the relativistic factor, which is employed during the particle pushing where v and c are the speeds of particle and light, u is the scalar momentum (u = γ v).
The planar finite strip-diode model is shown in Fig. 1. The 3-D sketch of the diode is exhibited in Fig. 1(a), where W is the width of both cathode and anode plates, D is the gap between the two plates, and U is the externally applied voltage. The electrons are emitted from the cathode more like a beam emitter without regard to the FE limit. As L (the length of the strip) is infinite along z axis, the 3-D model can be simplified to 2-D in Fig. 1(b). Therefore, the computational area is divided by uniform rectangles in 2-D Cartesian coordinate. The time step t obeys the Courant stability condition, where x and y are the lengths of the rectangle along x axis and y axis. As displayed in Fig. 1(b), in order to close the simulation area and ensure the convergence as well as the accuracy for solving Poisson's equation, the total computational space is set three times larger than the diode gap along both axes. The origin point is set at the lower left corner. To solve Poisson's equation, we should choose Dirichlet's or Neumann's boundary condition on the out boundary. For convenience, the homogeneous Neumann's boundary condition Eq. (6) is applied. Then, particle boundary condition and field boundary condition should both be considered in the simulation model. In most cases, the cathode and anode are made of metal. Hence, for the particle boundary, the primary electrons are all absorbed and the secondary electron emission could be ignored at the time the former ones arrive at the anode plate or return back to the cathode plate. What's more, the particles will be deleted immediately moving at the out boundary.
The fields at the two plates need to obey the metal boundary condition. For the field boundary condition at the out boundary, following Eq. (6) we can get the normal electric field condition. From the practical field distribution, the tangential fields are set the same at the upper and lower boundaries along y axis as Eq. (7), and are set the opposite at the left and right boundaries as Eq. (8).
∂ (x, y) ∂x ∂ (x, y) ∂y where (x, y)| y=3W and (x, y)| y=0 are at the upper and lower boundaries, and (x, y)| x=0 and (x, y)| x=3D are at the left and right boundaries, separately. Next, we apply the full electromagnetic PIC code UNIPIC [38]- [40], which is widely used for designing the vacuum electronic devices (such as high-power microwave [41], [42], millimeter wave [43] and terahertz devices [44], [45]), to study the SCL current density affected by the geometrical parameters gap separation D and width W , as well as the external parameters such as applied voltage U , the suppression magnetic field B, the initial energy distributions and temporal profiles of emission electrons.

III. EFFECTS OF GEOMETRICAL FACTOR, APPLIED VOLTAGE AND GUIDING MAGNETIC FIELD
A. THE OVER-SCL EFFECT Fig. 2 shows the currents on the cathode and anode plates respectively when the electrons are injected with a negligible initial energy. As we can see, the cathode current (red line) and the anode current (black line) will be of the same limited value at the late steady state, where all the average current densities in the diode gap are stable at the YOZ plane. We know from Ref. [11] that when the electron reversal occurs, the virtual cathode firstly forms and the increase rate of transported current density with the injection current will start to decline. When the electron reversal is as strong as the electron injection, the increase rate of current density declines to zero and the transported charge or current density reaches the maximal. For convenience, at the steady state, we choose the anode current density or cathode current density as the SCL current density or the largest average current density [21]. The anode current is acquired by computing the charge of electrons absorbed by anode per t, while the cathode current is obtained by computing the charge passing through the plane one grid away from the cathode along x axis per t. Also, it is observed that when the injection current is larger, the average current density is even going to decrease under the stronger impact of the self-induced electric field which pulls the electrons back, as compared between the blue line and the black line shown in Fig. 2, and we term it as the over-SCL effect, which is also demonstrated by Zhu et al. [46].  Fig. 3 not only proves to accord with the law in Ref. [21] and our previous work [18] that J CL (2)/J CL (1) grows larger with the decrease of W /D, but also generalizes the scale law to different gap separations. The three hollow symbols in Fig. 3 are the data in Ref. [21], which matches with our results well except for the extremely small W /D. It is obvious in Fig. 4 that since the Coulomb force from the particles are going to repel one another, when W /D is relatively small and suppression magnetic B is not strong enough, the moving particles will escape from the diode easily and periodically. As B is getting smaller, the escape spreads along x axis further. We can also find that when there is no suppression magnetic field, particles seem to distribute uniformly in the computation space. Except the suppression magnetic field, the transported current will generate selfinduced magnetic field which is known proportion to the radius of the electron beam. When W /D and W are small, the radius of the transported beam is small as well, and the selfinduced magnetic field is not strong, which makes it more likely to observe the electron escape. Additionally, the edgeeffects are quite significant in Figs. 4(a) and 4(b) that the majority of the particles gather around the boundary of the particle beam, so the density of transported electrons is low and even discrete around the central region of the diode. And the reason why the current in the center region is discrete is due to the periodical formation of virtual cathode, which limits the electron emission ability. Fig. 5 gives the particle distribution when the geometrical factor comparably large. It is observed that at this time a suppression magnetic field of 5 T will constrain most of the particles, and only when B approaches 0, will particles escape from the diode gradually.  Meanwhile, because the weak magnetic field is not able to constrain the particles, we can find from Tab. 1 that the anode current descends with the decrease of B. This discovery is opposite to Ref. [21] that the external magnetic bring little influence on the current. We guess the computation model in their simulation is relatively large (D = 0.01 m) so that negligible number of particles escape from the diode.  voltage. The four hollow symbols are the results in Ref. [21], which also match well with our computation results.   7 displays the normalized average current density J x , electric field E x , and E y along y axis in the diode. The J x and VOLUME 10, 2022 E x are symmetric about the central axis of the diode, and the maximum locates at the two edges, which is consistent with that in [24]. The distribution of the electron beam should be considered the competition of the Coulomb force and selfinduced or suppression magnetic field force. The electrons at the edge will be subjected to larger repulsive Coulomb forces from the inner electrons, and since the magnetic field will constrain the out-escaping electrons, the majority of electrons will converge around the edge, behaving as the edge-effect. As a result, larger electron density causes the larger electric fields at the edge. While, though the maximal y-direction electric field is distributed at the edge, the polarity is opposite at each edge, and that's the reason we use the symmetric boundary condition.  Fig. 8(a), we can find that when W /D is small and the suppression magnetic field is weak (B = 0.1 T) to restrain the particles, increasing the voltage will make the injection particles more energetic and spread more uniformly in the computation space. On the other hand, when the magnetic field is strong enough, the applied voltage will not impact the particle distribution but only increase the transported current, as in Figs. 8(b) and 8(c). What's more, comparing Figs. 8(b) and 8(c), it is indicated that the nonuniform and periodically discrete particle distribution occurs in smaller W /D condition because of the formation of virtual cathode.

C. EFFECT OF BIASED VOLTAGE UNDER VARYING SUPPRESSION MAGNETIC FIELD
In conclusion, the geometrical factor W /D, gap width D and suppression magnetic field B will impact the particle escape phenomenon. When they are all small, particles are observed to escape more easily because the constrain of magnetic field force on the particles cannot offset the Coulomb force from the inside particles. Meanwhile, significant edge-effect is exhibited that the majority of particles gather around the two edges of the current beam and thus form a thin and periodically discrete 'road' through the center region. The edge-effect is obvious after the formation of virtual cathode when W /D is relatively small but the magnetic field is strong enough. Two reasons accounting for the edge-effect: 1) The virtual cathode forms preferentially in the center region, which limits the electric field and electron injection greatly [10]; 2) The edge-effect enlarges the local electric field and electron injection intensity. Meanwhile, our work provides supplementary instruction about J CL (2)/J CL (1) varying with gap separation D and applied voltage U under strong magnetic field to Ref. [21]. Last but not the least, the particles escape and then the resulting less anode current are found in our work.

IV. COMBINED EFFECT OF ELECTRON EMISSION TEMPORAL PROFILE AND INITIAL VELOCITY
The CL law is used to describe the maximal current density at steady state only when the injection current is timeindependent. Nevertheless, if the injection current intensity changes with time, the SCL effect will change as well, as the SCL current density is influenced by the change rate of the injection current, i.e., by both the temporal profile and the injection current peak. Therefore, the way by only increasing the injection current peak in a time-dependent profile to find the maximal transported current density is incorrect.
In this case, Ref. [27] defined the average current density at the time the backflow firstly observed as the SCL current density of different temporal profiles. But in fact, it is quite a hard work to find the critical point. As shown in Fig. 2, larger injection current will decrease the transported current at some inflection point and the point is easier to obtain [46]. Thus, in this article, we choose the maximal transported current density as the value of interest to study. But due to the pulse width is not so long to obtain the steady state in diode, now, let us choose where we should choose the diagnostic point.
Tab. 2 exhibits the anode and cathode currents in the diode when the temporal profile is linearly increasing with the rise time (from zero to the peak) t r = 5 ns and t r = 0.005 ns which are larger and smaller than the transit time T CL [11] respectively.
where D = 9.9 × 10 −3 m, W = 2 × 10 −2 m (W /D = 2.02), L = 1 × 10 −3 m, U = 10 4 V and thus T CL = 0.5 ns. Selfinduced electric field in the diode will modulate the moving electrons especially for the short pulse beam, widening the pulse and shortening the peak of current from cathode to anode. Thus, the anode and cathode currents are not the same.

A. DEFINITION OF J max
For the modulation of electric field and the infinitely increasing currents with the injection current, we choose the average current density over cathode as the investigative maximal transported current density J max and set an injection current of a fixed peak I p = 5 A mainly in order to obtain the tendency of the change of J max under different temporal profiles combined with electron emission spectra. Two temporal profiles are set to compare the influence of linearly and non-linearly increasing, where I p = 5 A is the peak current. And Eq. (10.1) represents linear temporal profile [27], while Eq. (10.2) is the sinesquare temporal profile which is a classical temporal profile for photoelectron emission [47]. Fig. 9 shows the ratio of the maximal transported current density on the cathode J max to 2-D SCL current density J CL (2) varies with t r /T CL , where the temporal profiles are sinesquare and linear, respectively. It is indicated that there are three stages for the change of the current density, where the first two stages are similar to that in [11]. However, since the injection current is fixed 5 A rather than an adjustable value, the specific results cannot be the same as that in Ref. [11], the more we care for is the general tendency. For stage I, as we can see, when t r is larger than T CL , J max /J CL (2) will still change a little. It is because T CL is calculated with a timeindependent profile of which injection current is always the largest, while the transported charge for the time-dependent profile is smaller in one transit time starting from zero. Thus, the current density cannot be steady so soon and the SCL force is not so strong as the time-independent one. It is interesting that when t r decreases to some threshold in stage III, pointed by the rectangle in Fig. 9, the transported current will not reach the injection current peak as supposed but even decrease. At this time, the transported current cannot reach steady so that the current density will always change through the diode. Such few charges under the small t r cannot support the transport current to increase further. Though the injection current peak limits the maximal transported current density, we still need to fix the injection current because increasing the injection current peak will result in the same current increase rate by shortening the rise time.

C. EFFECTS OF INITIAL ENERGY SPECTRUM AND TEMPORAL PROFILE
The injection electrons with initial energy are more able to overcome the SCL force, and our former work has found the significant impact of the initial energy of different spectra on the SCL current density [18]. Now, we use the 2-D model with D = 9.9 × 10 −3 m, W = 2 × 10 −2 m (W /D = 2.02), L = 1 × 10 −3 m, to study the combined effect of temporal profile and initial energy. Average initial electron energy E av under a spectrum f (t) is described as We simulate the transported current peak varying with initial electron energy of three spectra, mono-energy (simplified model), exponential distribution (for photoelectron) and Boltzmann distribution (for thermal electron), expressed as Eq. (12).
where, f A in Eq. (12.b) is a normalization factor, k is the Boltzmann constant, E in Eq. (12.b) and (12.c) is from zero to infinity, and T is derived from E av = 1 2 mv 2 = 3 2 kT . Four temporal profiles are applied, namely sine-square and linear temporal profile with t r = 0.1 ns and t r = 0.5 ns separately. The injection current peak is chosen I p = 5 A. From Fig. 10 we can find that the spectra truly affect the behaviors of current density. Although the spectra and temporal profiles are different from one another, the general tendency of J max /J CL (2) vs. E av is similar. Since the proportion of high energy electrons changes a lot with E av under the three energy spectra, the ability to weaken SCL effect is different, and thus there exist some crosses between the curves. Then, some conclusions are listed as follows: (A) With the same temporal profile and initial electron energy E av , the shorter rise time gives the larger current density J max . (B) With the same t r under different temporal profiles, the nonlinear change of injection electron is going to bring the larger average current density J max . (C) With the same temporal profile, J max is more sensitive to E av when the rise time is larger. (D) Exponential distribution brings the less energies to injection electrons than others, and Boltzmann distribution gives larger J max than mono-energy at least when E av /eU < 10 −3 . Therefore, though we know the steady current density which transports electric quantity stably is limited by the CL law, if we only concentrate on the transient output, the unstable current density enhancement is acceptable as well. In order to enhance the SCL current, enlarging the increase rate of injection current with a proper rise time is a worthy choice, and a better electron injection mechanism should be considered. As it is tough to obtain the J max which varies with both temporal profile and initial energy of spectra, only some phenomenological results are discovered, and further researches need to be addressed.

V. CONCLUSION
This paper presents the SCL current density and particle distribution simulated by the UNIPIC code in a 2-D planar diode over a finite strip based on the Lau's model [19]. As a result of the significant edge-effect, larger current density and electric field are generated at the edge, leaving a thin and discrete (due to the Coulomb blockade) passage in the center. Moreover, considering the effect of geometrical parameters, periodical particle escaping from the diode gap is found when suppression magnetic field B, gap separation D and geometrical factor W /D are all small at the same time, leading to the decline of anode current density. The results also demonstrate that when B is strong enough to overcome the Coulomb force and restrain the particles in the diode, J CL (2)/J CL (1) is not sensitive to gap separation D and applied voltage U , which proves the law in Ref. [21] is of good applicability to different inputs.
As for time-independent injection, due to the SCL effect, the transported current will be modulated and the anode current density will be smaller than the injection one. Hence, increasing the injection current is an adoptable way to obtain a larger transported current. However, it is not a generally applied method, because the subsequently enormous SCL effect will decrease the transported current density. When it comes to a time-dependent injection, following Valfells's work, we find that increasing the change rate of the injection current by employing the nonlinear profile or decreasing the rise time no less than 10 −3 · T CL is the considerable method to obtain a larger transient current density. However, obeying the CL law, the stable transported current density is still limited. Meanwhile, the study on cylindrical or spherical models has been reported before [48], which is more valued in engineering application and deserves to be studied deeply. For an actual condition, thermionic emission could cause larger current density than photoemission under the same average electron energy. Since the actual electron injection is timedependent and energy-dependent, these parts deserve further studies. Last but not the least, as the field emission is a more common electron emission method, and the transition from FN law to CL law is practical [49], [50] especially under the rough surface condition [51], it is urged to take all the demands above into our codes in the future. He is currently working with NINT as an Associate Professor. His research interests include electromagnetic pulse simulation and plasma physics. He is currently working at NINT as a Research Assistant. His research interest includes numerical electromagnetic methods.