Relativistic time-dependent quantum dynamics across supercritical barriers for Klein-Gordon and Dirac particles

We investigate wavepacket dynamics across supercritical barriers for the Klein-Gordon and Dirac equations. Our treatment is based on a multiple scattering expansion (MSE). For spin-0 particles, the MSE diverges, rendering invalid the use of the usual connection formulas for the scattering basis functions. In a time-dependent formulation, the divergent character of the MSE naturally accounts for charge creation at the barrier boundaries. In the Dirac case, the MSE converges and no charge is created. We show that this time-dependent charge behavior dynamics can adequately explain the Klein paradox in a first quantized setting. We further compare our semi-analytical wavepacket approach to exact finite-difference solutions of the relativistic wave equations.


I. INTRODUCTION
One of the salient features of relativistic quantum mechanics is the intrinsic mixture of particles with their corresponding antiparticles. This aspect already appears in "first quantized" single particle Relativistic Quantum Mechanics (RQM) [1]. This phenomenon becomes prominent when a particle is placed in a very strong external classical field -a field strong enough to close the gap between particles and antiparticles, known as a supercritical field.
Most textbooks, as well as the overwhelming majority of past works employ timeindependent stationary phase or plane wave arguments when dealing with supercritical fields within RQM. As is well known, in quantum scattering theory time-independent quantities such as cross-sections or energy levels can be readily computed, but it is often ambiguous to attempt to infer the dynamics from considerations involving a single plane wave. This is even more the case when relativistic phenomena are investigated: one has to deal with additional difficulties, such as the breakdown of the usual connection formulas linking the scattering solutions in different regions, or the Klein paradox -a phenomenon usually expressed as a current reflected from a supercritical step or barrier higher than the incident one. Standard textbooks (e.g. [1][2][3][4]) give different, often conflicting accounts of the Klein paradox. This situation is reflected even in recent works [5][6][7][8][9][10][11][12][13][14][15], that reach different conclusions generally based on time independent considerations.
In part for this reason, it is often stated, in different situations dealing with the Klein paradox, that a RQM approach is inadequate and that a QFT treatment is necessary (see e.g. [3] for the step case, or [8] for the barrier case). It should be nevertheless recalled that the relativistic quantum mechanical equations have positive and negative energy states built-in. This is why although formally remaining a one-particle theory, this particle can be in a superposition of states with different charges (related to the particle and its associated anti-particle). The total charge will be conserved, but locally (in the cases we will be dealing with, on either side of regions in which the potential changes) positive or negative charges can increase locally and propagate. We will assess whether this leads to a consistent timedependent first quantized explanation of the Klein paradox. Interestingly, other very recent works have focused on different time-dependent aspects of the relativistic wave equations [16,17].
In this paper we will develop a time-dependent wavepacket treatment suited to investigate the dynamics of spin-0 bosons and spin-1/2 fermions in model supercriticial barriers. In order to develop our semi-analytic wavepacket approach, we will rely on a Multiple Scattering Expansion (MSE). We will see that the nature of the MSE is different for solutions of the Klein-Gordon and Dirac equations. In the Klein-Gordon case, the MSE diverges, physically corresponding to charge creation as the wavepacket hits a barrier's edge. This implies that the usual connection formulas between wavefunctions in different regions -which are obeyed when a single step potential is considered -are not applicable when dealing with supercritical barriers. Employing connection formulas leads to inconsistencies, like superluminal traversal times, that were recently noted [8] though (incorrectly, as we will see) attributed to the limitations of the first quantized formalism. This is why it is crucial to take into account the multiple scattering process inside the barrier. In the Dirac case, the MSE converges, as no charge is created when the wavepacket hits the barrier. This leads to an interpretation of Klein tunneling that is qualitatively different from the bosonic case. This also ensures that the connection formulas, that have been employed in a countless number of works, remain valid. Our semi-analytical wavepacket treatment will be compared for both the Klein-Gordon and Dirac cases to exact numerical results obtained from finite-difference computations.
The multiple scattering expansion will be given in Sec. II, after introducing the model barriers we will be working with. This will be employed in Sec. III to build wavepackets. We will first detail the construction of wavepackets from plane-wave expansions weighted by the scattering amplitudes. We will then give a couple of examples displaying the dynamics of a bosonic or fermionic particle impinging on a supercritical barrier. Our wavepacket results will be compared to numerical solutions obtained from a code we have developed to solve numerically the relativistic wave equations. Our results will be discussed in Sec. IV. We will more specifically focus on the extent to which the Klein paradox can be accounted for within a first quantized framework. We close the paper by a Conclusion (Sec. V), and two Appendices give additional detail on the multiple scattering expansion and on the finitedifference scheme employed for the numerical computations.

A. Relativistic wave equations
For definitiness and future reference, let us recall the wave equations for a "particle" impinging on a one dimensional barrier. Note we will mostly be dealing in this paper with a barrier, that we discriminate from a step, which is by far the case that has most often been considered in studies of the Klein paradox. In the barrier case, the potential is appreciably different from zero in a region 0 < x < L, where L is the barrier width. The simplest case to treat is the rectangular barrier defined by where θ denotes the unit-step function and V denotes the barrier height. We will also consider smooth barriers, for which we will employ the potential since for this potential analytical solutions are known [18]; is the smoothness parameter.
The Klein-Gordon (KG) equation, describing spin-0 bosons of rest mass m, in the presence of a time-independent potential (electrostatic type) in the minimal coupling scheme is given in the canonical form by We will mostly use the Hamiltonian version of the KG equation, given by where Ψ has components linked to the solutions ψ(t, x) of the canonical KG equation through σ j with j = 1, 2, 3 denotes the Pauli matrices. Recall [1] that the local charge probability density ρ(t, x) can be negative; it is given in the Hamiltonian formulation by so that the positive and negative charge amplitudes are respectively carried by ϕ and χ.
The more general expression for the Klein-Gordon scalar product for two wavefunctions Ψ 1 = (ϕ 1 , χ 1 ) and Ψ 2 = (ϕ 2 , χ 2 ) given by Eq. (5) is given by The Dirac equation for a fermion propagating in one spatial direction can be reduced to an equation for a 2 component wavefunction Φ provided the spin is "frozen" (with no spin degree of freedom, in particular no spin flip at the barrier). It can be seen either by starting from the usual Dirac equation with an electrostatic potential, where γ µ , µ = 0, ..., 3 are the usual 4 × 4 gamma matrices [1], and then restricting the 4 component wavefunction Υ(t, x) to the non-trivial subspace involving only two components (valid in the laboratory frame); alternatively one can start from a non-covariant equation implementing the usual constraints leading to the Dirac equation for one spatial dimension [15]. This yields the 2-component Dirac equation for Φ Recall that in the Dirac case the local probability density ρ(t, x) = Φ † (t, x)Φ(t, x) is always positive and can be written as ρ(t, x) = |ϕ(t, x)| 2 + |χ(t, x)| 2 if we label the components as in Eq. (5). The positive definite scalar product is

B. Supercritical barriers
A supercritical potential is a potential high enough to give rise to Klein tunneling [19], whereby the incoming wavepacket penetrates undamped inside the barrier. This happens when the potential is high enough (and its slope is steep enough) to bridge the gap between positive and negative energy solutions. In the bosonic case, this gives rise [20] to superradiance (a reflected current higher than the incoming one). In the fermionic case there is no superradiance (although some authors suggest differently, see e.g. [2,13]), and supercritical steps have been deemed to have an acceptable interpretation only within a Quantum Field Theory (QFT) approach [21,22], a point that appears to be supported by the wide variety of conflicting interpretations of Klein tunneling that have been proposed within the first quantized framework [5][6][7][8][9][10][11][12][13][14][15].
As is well known [1,20], relativistic scattering of plane waves on a subcritical potential is similar to the usual non-relativistic situation: for a sufficiently high barrier, the transmission amplitude is small and tunneling takes place through exponentially decreasing waves. For example for the rectangular potential given by Eq. (1), plane wave solutions of the canonical KG equation are of the form where j = 1, 2, 3 denotes the regions depicted in Fig. 1 and the ± signs corresponds to states with positive and negative energies respectively. p j is the momentum; for positive energies, p j > 0 gives a wave moving from left to right (but from right to left for negative energies).
Assume boundary conditions for which an incident positive energy wave travels from left to right; this imposes B 3 = 0 and for definiteness we set the incoming amplitude to A 1 = 1.
The other amplitudes A j and B j are deduced by matching the wavefunctions and their space derivatives at the boundaries x = 0 and x = L (for reasons that will become clear below, we will not need to deal with boundary conditions for negative energy plane waves; we henceforth write A for A + , etc.). In particular the reflection amplitude r = B 1 and the transmission amplitude t = A 3 are obtained by employing these connection formulas matching the wavefunction and the first derivative at the boundaries (in the Dirac case, only the two components of the wavefunction are matched, since the equation involves only a first space derivative). The connection formulas do not necessarily hold when V becomes supercritical, as we will now discuss.

C. Klein-Gordon equation: Divergent Multiple Scattering Expansion
It is well-known that when transmission of waves across several media takes place, one has to take into account a multiple scattering process. Referring again to Fig. 1, a wave incoming from the left towards the barrier first encounters the left edge of the barrier. We set again A 1 = 1, and complete the boundary conditions by requiring that there is no wave coming towards the barrier from the right (so that at x = +∞ there is only a transmitted wave), setting B 3 = 0.
Consider therefore the asymptotically free (at x = −∞) wave coming toward the barrier.
Reflection takes place with amplitude r i l , which is the reflection amplitude of a step. The transmitted amplitude at that point, t i l is that of a step, but the wave traveling inside the barrier reaches the right edge, at which part of it gets transmitted, with amplitude t i r , and part reflected with amplitude r i r . This reflected wave travels back towards the left side of the barrier, getting reflected and transmitted with coefficients r o l and t o l . This process is iterated an infinite number of times, and the overall reflection and transmission amplitudes, B 2 and A 3 (as well as the amplitudes inside the barrier, at least in the rectangular barrier case) appear as The amplitudes obtained by using this Multiple Scattering Expansion (MSE, also called Multiple Reflection Expansion [23] ) approach should match those obtained by employing the usual connection formulas at the boundaries, but this will happen only provided the sums in Eq. (15) converge.
Let us consider the Klein-Gordon equation in the presence of a supercritical step potential at x = 0, Scattering of a positive energy plane wave coming from the left, as given by Eq. (13) sets since there is no multiple scattering for a simple potential step (we use the bar to avoid confusion between the step and barrier amplitudes).Ā 2 andB 1 are thus obtained by applying the boundary conditions at x = 0, giving the well-known result where following our notation given above, p 1 > 0 is the momentum of the incoming plane wave (in region 1) and Note that p 2 must be negative in order to fulfill the boundary conditions, since a negative energy plane wave with positive momentum moves from right to left [1,20].
The amplitudesB 1 andĀ 2 of the step correspond to the amplitudes r i l and t i l of the barrier MSE, Eq. (15). The other amplitudes are obtained by considering a step described by the potential V θ (L − x) (this gives r i r and t i r ) or the same step (16) but with a wave coming from the right ("outgoing" boundary condition) at x = 0 (to get r o l and t o l ), as detailed in Appendix A. By inserting the values for these elementary scattering amplitudes into the MSE of Eq. (15), we obtain the scattering amplitudes for the rectangular barrier. These in which case it can be checked, using (1 − ζ) −1 = ∞ x=0 ζ n , that the amplitudes A j and B j match the ones obtained by using the connection formulas at the boundaries.
The condition (19) is fulfilled for the KG equation for subcritical barriers, since then both p 1 and p 2 are positive. However for supercritical barriers, we have p 1 > 0 and p 2 < 0 and the MSE diverges. This divergence doesn't make much sense in a stationary plane-wave picture, in which the scattering amplitudes become infinite, but we will see below in Sec. III that in a time-dependent approach, the divergence corresponds physically to the creation of charge each time a wavepacket hits a barrier edge. The negative charge inside the barrier increases, and a corresponding positive charge is created outside the barrier, so that the total charge is conserved.
Interestingly, if the "converged" amplitudes (usually obtained by employing the connection formulas) are employed in the supercritical case, unphysical results are obtained. It was recently noticed [8] in a plane wave analysis that the barrier traversal time defined from the phase energy derivative was superluminal in the supercritical case, an unphysical result attributed to the limitations of the first quantized formalism. In a wavepacket approach, building wavepackets with the converged amplitudes results in an acausal wavepacket coming out from the right side of the supercritical barrier before the incoming wavepacket has even hit the barrier [23]. Other works have also employed the connection formulas in a supercritical context [6,12].

D. Dirac equation: Convergent Multiple Scattering Expansion
The derivation of the MSE for the Dirac equation (11) is similar to the one we have given for the KG equation. The plane wave solutions of the 2-component state Φ(t, x) take the form: where j refers again to the 3 regions depicted in Fig. 1 and n ± j are normalization coefficients. As in the KG case, let us treat the barrier as a multiple scattering expansion, with the same boundary conditions The amplitudes A ± j and B ± j are again given by Eq. (15), with coefficients r i l , t i l ... that are different for positive and negative energy wavefunctions. As in the KG case, the MSE for a rectangular barrier is built from the reflection and transmission amplitudesĀ ± j andB ± j of the steps V θ(x) and V θ(L − x). Since the Dirac equation (11) is first order in the space variable, the positive energy amplitudes for j = 1, 2 are obtained from the matching condition Φ + 1 (t, The amplitudesB 1 andĀ 2 of the step correspond to the amplitudes r i+ l and t i+ l entering the MSE for the barrier. The other step amplitudes are obtained in the same manner; they are given in Appendix A. The barrier amplitudes A ± j and B ± j can now be obtained from the MSE (15). It can be seen that the series converge provided which is verified because α + 2 (p 1 ) is positive since V is supercritical. Therefore for the Dirac equation, the MSE converges. The usual connection formulas may be employed to obtain directly the barrier amplitudes A ± j and B ± j (for the rectangular barrier, these are Φ ± 1 (t, ; in the general case, the solutions in the barrier region are matched to the asymptotic scattering solutions). Note that most past works (e.g., [12,24,25]) have indeed employed the connection formulas in this manner, without however examining the justifications for their use. Note also that this implies the different behavior fermions have relative to bosons when scattering on a supercritical barrier: superradiance only happens for bosons, and this is due to the divergent character of the multiple scattering expansion (15).

A. Construction from plane-wave expansions
The most straightforward way to construct wavepackets starting from an initial distribution is to employ a plane-wave expansion, valid asymptotically far enough from the barrier.
For a sufficiently steep barrier with a nearly constant value (see Fig. 1), the plane-wave expansion is valid everywhere except in the slope region. For a rectangular barrier, a planewave expansion is exact everywhere.

Klein-Gordon equation wavepackets
The Klein-Gordon plane-waves were given in Eq. (13). These solutions can be expressed in the Hamilton form, Eq. (5), as: where p 3 = p 1 and p 2 is given by Eq. (18), the amplitudes A ± j and B ± j are given by Eq. (15), and N is a global normalisation constant. To simplify the notation, we will put and we define B ± j similarly. Assume that at t = 0 we have an initial wavefunction G(t = 0, x) = (ϕ G (0, x), χ G (0, x)) in region 1, to the left of the barrier. The time evolution can be employed by applying the pseudo-unitary evolution operator on G(t = 0, x), or equivalently by using the Fourier transform G(t = 0, x) = dpe ipx/ Ĝ (t = 0, p). In both cases we need the projections c ± KG = Ψ ± 1 G KG on the plane wave solutions Ψ ± 1 (0, x) given by Eq. (27). The c ± KG are obtained by using Eq. (9), and the time evolved wavepacket can then be written as where θ j ensures each expression is used only in the region j in which it is valid, as per Fig.   1 (explicitly, θ 1 = θ(−x), θ 2 = θ(x)θ(L − x) and θ 3 = θ(x − L)).
To be specific, we will choose an initially localized state of the form (we have omitted an overall constant term depending on the normalization and Gaussian width), we see that c − KG is non-zero (although it is small for non-relativistic velocities). This means that there is a negative energy component that develops for t > 0. It turns out that this negative energy component moves to the left (recall that p 0 > 0 yields an antiparticle that moves in the negative direction), so only the positive energy wavepacket (particle) impinges on the barrier. We might as well have chosen to restrict G(0, x) to the expansion over the positive energy eigenstates in Eq. (29), but when comparing with numerical solutions it is more interesting to keep (30) as the initial state. The wavepacket is therefore given by where the amplitudes A j and B j are given by Eqs. (15) from the coefficients collected in Appendix A, and the unit-step function are defined below Eq. (29).

Dirac equation wavepackets
A similar construction can be used to build wavepackets evolving according to the Dirac equation, starting from an initial state |G(t = 0) that is expanded over the Dirac equation plane-wave solutions given by Eq. (20). The coefficients c ± D = Φ ± 1 G are now obtained by using the standard positive definite scalar product (12). For definiteness, let us start again with an initial state G(t = 0, x) of the same form as Eq. (30). This wavepacket will have projections on the plane wave solutions with positive and negative energies (the latter vanishing only in the non-relativistic limit). The coefficients are given by and we set n ± (p 1 ) = 1 + α ± 1 (p 1 ) 2 −1/2 . The time evolved wavepacket written in a compact form is Since the MSE converges, the amplitudes A ± j and B ± j are most easily obtained by the matching conditions as explained below Eq. (26) (see Appendix A). As in the Klein-Gordon case only the positive energy contributions impinge on the barrier, as the negative energy contributions move to the left.

B. Numerical solutions
We have computed numerical solutions to the relativistic wave equations (4) and (11). This was done by discretizing the corresponding evolution operator in real space for small time steps. The initial wavefunction G(t = 0, x) is discretized on a fixed space-grid and the derivatives in the evolution operator are approximated by finite differences in the fourth or fifth-order approximation (see Appendix B for additional details).
For the Klein-Gordon equation, real-space approaches are now possible [26,27] although they are computationally more demanding than the Fourier-transformed based split operator methods. The evolution operator for an infinitesimal time step δt is given by is obtained from the right hand-side of Eq. (4). This operator is pseudo-unitary in the sense that σ 3 U † KG σ 3 U KG = I. It is given explicitly by where u KG is given by Eq. (B3) and n max sets the number of terms kept in the Taylor expansion.
Similarly, an initial Dirac wavepacket G(t = 0, x) is evolved numerically according to Eq.
(11) by implementing the finite-difference approximation of the Hamiltonian The time evolution operator is then obtained as where u D is given by Eq. (B8). Note that the 'fermion-doubling" problem [28] (appearance of spurious numerical solutions to the Dirac equation) is avoided by working with small space steps so that the momentum cutoff remains well below the largest momentum included in the wavepackets. Further computational details are collected in Appendix B.

C. Illustrative Results
We will now illustrate our wavepacket approach to relativistic scattering on a barrier, and compare it with fully numerical solutions obtained by solving the relativistic wave equations with a finite-difference scheme. An interesting case to consider is a supercritical barrier (since we can test the divergent and convergent characters of the Klein-Gordon and Dirac equations respectively), with a flat summit (so that the plane-wave expansion remains valid and can be tested against the numerical computations inside the barrier). A rectangular barrier, or a smooth version, given by Eq. (2) above, is well suited to our investigation. We pick the initial state Eq. (30) with x 0 = −400, p 0 = √ 5/2 and d = 50 (in order to have a rather narrow momentum distribution). We also provide numerical solutions obtained In practice, the number of terms that need to be taken into account in the MSE sum (15) is congruent with the time t at which the wavepacket is computed. Indeed, the nth term corresponds formally to a wavepacket translated by the order of 2nπL, that for values of n that are high enough did not have time to reach the barrier. Hence these terms do not contribute to the wavepacket (in the calculations shown on Fig. 2, we included terms up to n = 4, since including additional terms did not make any difference). Note also that the particle and antiparticle character of the wavepacket can be inferred by looking at the component decomposition of G(t) = (ϕ G (t), χ G (t)) . Inside the barrier, we have G (0, χ G ) , i.e. mostly negative charge, while outside the barrier, G (ϕ G , 0) (mostly positive charge) as pictured in Fig. 3. The wavepacket dynamics for a spin 1/2 fermion is shown in Fig. 4. We picked the same barrier and initial state parameters as in the Klein-Gordon case except that the smoothness parameter was taken to be the size of the space integration step in the numerical code ( 250), effectively corresponding to the rectangular barrier limit. We have therefore employed the amplitudes (A9) of the rectangular case. We note again the very good agreement between the wavepackets constructed with Eq. (34) and the numerical solutions. Klein tunneling now takes place without charge creation, and the wavepacket inside the barrier disappears after a few oscillations.

A. General remarks
We have developed a wave-packet approach able to describe Klein tunneling across a supercritical barrier both for Klein-Gordon and Dirac particles. The main ingredient was seen to be the multiple scattering expansion, that diverges in the KG case but converges in the Dirac case. This gives a very different behavior for bosonic and fermionic wavepackets: the Klein-Gordon equation contains by construction charge creation. In the single particle framework this should be interpreted as a single particle being in a superposition of positive and negative energy eigenstates, though outside the barrier we essentially have a particle, while inside the barrier an anti-particle. The total charge is conserved, but inside the barrier the negative charge grows, while the positive charge increases outside the barrier (see Fig.   3). The Dirac equation on the other hand has no provision for charge creation. The charge of the incoming particle (say, electron) is distributed among the reflected, transmitted and inner barrier wavepackets.
The main issue now is whether these results lend to a consistent interpretation of the Klein paradox within the first quantized formalism. This is known to be problematic when basing considerations on stationary plane waves, and indeed different, sometimes conflicting interpretations of superradiance and supercritical tunneling in a first quantized framework have been given [5][6][7][8][9][10][11][12][13][14][15]. Klein tunneling has even been denied to exist in the bosonic case [14] or for fermions [10]. Of course since particle creation is induced by supercritical potentials, a proper approach requires a Quantum Field Theory (QFT) treatment. But time-dependent QFT approaches to tackle this problem are scarce (with the exception of the work reviewed in Ref [22], where supercriticial steps, rather than barriers, were investigated both for fermions and bosons; see also [29] and Refs. therein for recent ramifications of this time-dependent QFT method). Hence, although time-independent approaches to Klein tunneling within the first quantized framework might be seen as unreliable, it would be worthwhile to examine the situation with a time-dependent first quantized account, and examine its consistency with time-dependent QFT approaches.
B. Klein Paradox in a first quantized framework

Bosons
For the Klein-Gordon equation, the account seems rather simple: a supercritical potential creates positive and negative charges in equal amounts (since the total charge is conserved) but in different spatial regions (in either side of the region in which the potential changes).
In our wavepacket approach, this is explained by the divergent character of the multiple scattering expansion, although this property is encapsulated in the pseudo-unitary evolution operator that is solved numerically in the finite-difference code.
Hence, when a particle impinges on a supercritical barrier, the reflected wavepacket has a higher amplitude than the incoming one (which is the signature of superradiance). The positive charge created outside the left edge of the barrier is compensated by the creation of negative charge inside the barrier (see Fig. 3), giving rise to an antiparticle wavepacket inside the barrier (Klein tunneling). A similar process takes place when the antiparticle wavepacket reaches the right edge of the barrier: a positively charged particle wavepacket leaves the barrier, while this is compensated by the creation of negative charge inside the barrier, so that the reflected antiparticle wavepacket's amplitude has increased. This process, with positively charged wavepackets going out of the barrier with increasing amplitude while the negative charge inside the barrier oscillates with an increasing amplitude goes on indefinitely (to the extent that the supercritical potential can be maintained despite the growing charge inside the barrier).
It is difficult to accommodate the charge creation mechanism by supercritical fields with the idea that the first quantized framework would describe a single particle in a superposition of states of different charges, with the charge in distinct spatial regions increasing with time (rather than explaining this feature in terms of particle-antiparticle creation). However, leaving this important physical issue aside, we see that the time-dependent wavepackets by themselves account for the bosonic Klein paradox -superradiance and antiparticle tunneling.

Fermions
For the Dirac equation, the dynamics is very different. Assume the incoming particle is an electron (with negative charge). Part of the incoming wavepacket gets reflected by the barrier, and part gets transmitted. Contrary to the bosonic case, there is no charge creation (hence no superradiance), so the reflected wavepacket has a smaller amplitude than the incoming one. Since the total probability density is conserved, requiring charge conservation implies that the wavepacket tunneling inside the barrier also has negative charge. This is difficult to accommodate within a single particle picture, even by relying on hole theory. propagation of such an electron is difficult to account within hole theory, but we will see that the propagation of a negative charge inside the barrier fits well with the physics afforded by a QFT account of supercritical potentials, as we now discuss.

C. Relation to time-dependent QFT
Klein tunneling and the Klein paradox ultimately depend on a multiparticle process involving particle-antiparticle pairs creation, that should therefore be described by Quantum Field Theory (QFT). A time-dependent QFT scheme [22] has investigated Klein tunneling for rectangular and smooth steps. The QFT calculations have been compared to numerical solutions of the Klein-Gordon [30] and Dirac [31] equations.

Bosonic QFT
In the absence of any incoming particle, space-time dependent charge densities obtained from the bosonic field operator indicate that the supercritical step induces pair creation, particles with positive energies moving away from the step and wavepackets with negative energies (antiparticles) "tunneling" inside the step [30]. When a particle is sent towards the barrier, the QFT computations show that pair creation is enhanced (precisely by the transmission amplitude t i l of Eq. (15)). This enhancement corresponds to the first quantized computations for the step. The same correspondence between QFT and our results can be thought to hold in the barrier case: the wavepackets we obtained would correspond to the pair creation enhancement produced by sending a particle on the barrier, on top of the pair creation process out of the vacuum. Note that in the barrier case, spontaneous pair creations are also expected to lead to an amplification mechanism through multiple reflections inside the barrier.

Fermionic QFT
In the Dirac case, QFT computations of the time-dependent spatial densities from the fermionic operators for the step given by Eq. (16) show that an incoming electron modifies the pair production process induced by the supercritical field [31]. The reflected fraction of the incoming first quantized wavepacket appears as an excess of the particle charge (relative to the charge produced by the step). The transmitted wavepacket, propagating inside the step, appears instead as a dip in the anti-particle (positronic) charge produced by the supercritical potential. The interpretation is that, as a result of the Pauli principle, the incoming electron decreases the pair-production process that takes place at the step. The decrease of antiparticle production at the edge of the step gives rise to a hole in the positronic charge inside the step. The dip in the electron production is partially compensated by the incoming electron that is interpreted as being fully reflected [31], yielding overall an excess of electronic charge which is the reflected first quantized wavepacket.
It is not obvious whether one can extrapolate straightforwardly the QFT results for a step to a fermionic barrier. Now both edges of the barrier have a pair production process. Inside the barrier, the Pauli principle also applies to positrons. This has no counterpart in a first quantized framework. Nevertheless, one can speculate that, at least for a sufficiently wide barrier, an incoming electron partially suppresses pair production (similarly to the step) at

V. CONCLUSION
We have investigated in this paper the Klein paradox for barriers from a time-dependent perspective within the first quantized framework. We have developed a semi-analytical wavepacket approach relying on the properties of a multiple scattering process inside the barrier. This yields a very different behavior for bosons and fermions. In the bosonic case, each collision of the wavepacket on an edge of the supercritical potential creates charge (as the multiple scattering process diverges), leading to superradiance. In the Dirac case Klein tunneling occurs without superradiance (the MSE converges). The wavepacket calculations were complemented with exact numerical solutions obtained by implementing a finite-difference code, leading to an excellent agreement for rectangular and smooth barriers.
We have argued that while a stationary first quantized approach to the Klein paradox has resulted in different and conflicting interpretations, a time-dependent account describes the dynamics without any ambiguities. We further believe such wavepacket calculations might be valuable in order to have a qualitative or quantitative understanding for processes that should in principle be described by spacetime QFT approaches, which are computationally much more involved.
with µ, ν and λ given by

Dirac Equation
One can proceed similarly for the Dirac equation. For a rectangular barrier it suffices to match the wavefunction at x = 0 and x = L (since the Dirac equation is of first order in x) for the left step, and then the right step. This yields, for positive energy waves (we drop the "+" superscript) Negative energy amplitudes can be obtained similarly (though we do not need them in this work).

Appendix B: Finite-Difference Computations
To solve the Klein-Gordon equation, Eq. (4), numerically, we use the finite-difference fourth order approximation of the second derivative with respect to position. This approximation is given by: The one dimensional space of length l is discretized into a lattice of N x points, X = {− l 2 + nδx ; n = 0, 1, 2, . . . N x} with δx = l/N x . Thus, this second derivative is represented by an N x × N x matrix: The potential V (x) is represented by a diagonal matrix N x × N x where the diagonal V (X) is obtained by evaluating the potential V (x) at every point of the lattice X. We then use Eq. (B7) and V (X) to construct the operator Using u KG and the matrix: that commutes with both ∂ 2 x and V (X), we obtain the time evolution operator U KG (δt) = e −iĤ KG δt given by Eq. (37). This operator is represented by a 2N x × 2N x matrix that is applied recursively on the one-dimensional vector Ψ(t, X) to obtain Ψ(t + dt, X) following The numerical stability criteria was implemented similarly as in Ref. [27].
The numerical solutions to the Dirac equation were obtained along the same procedure.
The time evolution operator, Eq. (39) was built using the fifth-order finite-difference approximation of the first spatial derivative This approximation is employed on a discretized lattice of dimension l, yielding the matrix: Along with the potential V (X) (that is again diagonal) we obtain the matrix representing the operator u D , and finally the Dirac time evolution operator as per Eq. (39).
In the computations displayed on Fig. 2, we used a lattice of length l = 4000 and N x = 3 × 10 6 points. The time step was δt = 10 −3 . For the finite-difference solutions to the Dirac equation displayed in Fig. 4 we used l = 4000, N x = 5 × 10 5 and δt = 2 × 10 −2 . The momentum cutoff, π/(l/N x ), is 2 orders of magnitude larger than the highest momentum contributing to the wavepacket.