Transverse structure of electron in momentum space in basis light-front quantization

We investigate the leading-twist transverse momentum-dependent distribution functions (TMDs) for a physical electron, a spin-1/2 composite system consisting of a bare electron and a photon, using the Basis Light-front Quantization (BLFQ) framework. The light-front wave functions of the physical electron are obtained from the eigenvectors of the light-front QED Hamiltonian. We evaluate the TMDs using the overlaps of the light-front wave functions. The BLFQ results are found to be in excellent agreement with those TMDs calculated using lowest-order perturbation theory.


I. INTRODUCTION
A key tool for studying hadron structure is the deep inelastic scattering (DIS) process, where individual partons are resolved. One can extract the Parton Distribution Functions (PDFs) [1][2][3][4] from such experiments. PDFs encode the distribution of longitudinal momenta and polarizations carried by the partons. Being functions of longitudinal momentum fraction (x) only, they do not provide knowledge about the transverse spatial location and motion of the constituents inside the hadrons. On the other hand, the Transverse Momentumdependent Parton Distribution Functions (TMDs) [5][6][7] provide essential information about the distributions of both the longitudinal momentum fraction (x) and the relative transverse momentum (k ⊥ ) as well as the orbital motion of partons inside hadrons, allowing us to draw three-dimensional pictures of the hadrons. They appear in the description of semi-inclusive reactions like the semi-inclusive deep inelastic scattering (SIDIS) [8,9] and Drell-Yan process [10][11][12][13][14]. Quantum Chromodynamics (QCD) factorization theorems allow one to relate physical observables, such as cross sections, to TMDs via perturbatively calculable kernels. Many theoretical studies have concentrated on the formalism and the accompanying factorization theorem of the TMDs [15][16][17][18][19][20].
At leading-twist, there are eight TMDs for a composite spin-1/2 system. They are characterized by different spin-spin and spin-orbit correlations of the partons and the system. Three of them are generalizations of the three familiar collinear PDFs, namely the unpolarized f 1 (x), helicity g 1 (x), and the transversity h 1 (x) distributions. TMDs reduce to collinear PDFs by integrating over the transverse momenta. Basically, TMDs are the * huzhi@impcas.ac.cn † xsq234@impcas.ac.cn ‡ mondal@impcas.ac.cn § xbzhao@impcas.ac.cn ¶ jvary@iastate.edu extended version of collinear PDFs and one can extract the parton densities for various configurations of the target and the parton polarizations from TMDs. TMDs are very useful to explain a range of phenomena even when we only consider the leading-twist factorization. For example, one can explain the Collins asymmetry using transversity TMD h 1 (x, k ⊥ ) [21][22][23]; one can explain the double spin asymmetry A LT in SIDIS [24] using g 1T (x, k ⊥ ), which encodes the information of a longitudinally polarized parton in a transversely polarized hadron; one can employ the distribution corresponding to a transversely polarized parton in a transversely polarized hadron, h ⊥ 1T (x, k ⊥ ), also known as pretzelosity TMD, to explain the A sin(3φ h −φ S ) U T single spin asymmetry [25].
TMDs being non-perturbative in nature are always difficult to obtain from the first principles of QCD and thus, they have been investigated in several QCD inspired models. Examples include the diquark spectator model [26,27], MIT bag model [28][29][30], covariant parton model [31], light-cone constituent quark model [32] and light-front quark-diquark model with AdS/QCD predicted wave functions [33,34]. Meanwhile, promising theoretical frameworks for illuminating TMDs also include the Dyson-Schwinger equations approach [35,36] and the discretized space-time euclidean lattice [37][38][39][40][41]. However, both the Dyson-Schwinger equations and the lattice QCD approaches are mostly working in the euclidean spacetime and therefore encounter challenges in accessing TMDs directly.
In recent years, Basis Light-front Quantization (BLFQ) [42,43] has been developed as a nonperturbative approach which is very useful for studying the structure of relativistic bound states [44][45][46][47][48][49]. Starting with the light-front Hamiltonian obtained from the Lagrangian via the Legendre transformation, we solve for the light-front wave functions (LFWFs). Once we have the LFWFs, we are able to investigate several properties of the physical system, including the anomalous magnetic moment and GPDs of electrons [44,50], internal structure of baryons [48] and mesons [49,51], as well as scattering phenomena [52][53][54].
In Quantum Electrodynamics (QED), while the electron is treated as an elementary field, it can fluctuate into an electron-photon pair, e.g., e → eγ → e, with the same quantum number. The virtual photon can further break up into virtual electron and positron pairs with all possible combinations. Thus, an electron is no longer an isolated particle but is surrounded by virtual clouds of electrons, positrons and photons. In analogy to the partonic structure of hadrons, one can interpret them as partons contained in the physical electron. TMDs of the physical electron, as a system composed of a bare electron and a bare photon, have been investigated in Ref. [55] using the lowest-order perturbation theory. Many other observables of the electron, e.g., the electromagnetic and gravitational form factors as well as spin and orbital angular momentum [56], generalized parton distributions (GPDs) [57][58][59][60],Wigner distributions [61], etc. have been investigated in a model based on the quantum fluctuations of the electron in QED.
In this paper, we employ BLFQ to obtain TMDs of the physical spin-1/2 particle in QED, i.e. electron. Only the leading two Fock sectors (| e and | eγ ) are considered. During the diagonalization of the Hamiltonian, we use a mass counter-term in the first Fock sector to renormalize the mass of the system to that of the physical electron and we obtain the LFWFs of the ground state for the physical electron. It is reasonable to treat the perturbative results as equivalent to experimental data in QED. We then compare the BLFQ results with the perturbative results reported in Ref. [55]. It is our aim to demonstrate that the BLFQ approach, which is inherently non-perturbative, provides results, up to finite basis artifacts, that agree with the perturbative results. We also aim to show how we may mitigate the finite basis artifacts by achieving convergence of suitably smoothed BLFQ results. In achieving these aims, we establish the methods to be used in applications that require a nonperturbative approach.
For the purposes of this work, we set the gauge link generating from infinite many photon exchanges in the final-state interaction [62,63] to unity. This leaves us five non-zero TMDs out of the eight leading-twist TMDs. We also present spin-spin corrections between the bare electron and the physical electron for different polarization configurations. A detailed analysis of the TMDs with a nontrivial gauge link in the BLFQ approach will be reported in future studies.
The paper is organized as follows. We introduce the basics of the BLFQ approach in section II. A discussion of mass renormalization is also given in this section. In section III, the TMDs of electrons are evaluated. The numerical results of the TMDs and the spin-densities in the transverse momentum space for different spin configurations of the bare and the physical electrons are presented and discussed in section IV. The summary is given in section V.

II. BASIS LIGHT-FRONT QUANTIZATION
BLFQ is based on the Hamiltonian formalism and adopts light-front quantization [64]. The mass spectrum and light-front Fock state wave functions are obtained from the solution of the following eigenvalue equation: where the Fock space representation of the system state with momentum P and light-front helicity Λ, |P, Λ is expanded in terms of multi-particle light-front basis states as [7], Here x i = k + i /P + is the longitudinal fraction, k ⊥ i is the relative transverse momentum and λ i is the lightfront helicity, both of the i-th constituent. n represents the number of particles in a Fock state. The physical transverse momenta of the constituents are p ⊥ i = x i P ⊥ + k ⊥ i and the physical longitudinal momenta of the constituents p + i = k + i = x i P + . The boost invariant lightfront wave functions ψ Λ λ1...λn are independent of the total momentum of the state, P + and P ⊥ and depend only on x i and k ⊥ i . Like in QCD, we treat the physical electron as a composite particle, where we could find the bare electrons, positrons and photons as its partons. In this work, we restrict ourselves to only the first two Fock sectors: The full Hamiltonian we diagonalize iŝ where the light-front 1 QED HamiltonianĤ QED = P +P − QED − (P ⊥ ) 2 . We obtainP − QED from the QED Lagrangian [64,65] through the Legendre transform and adopting the light-cone gauge A + = 0. Restricted to the leading two Fock sectors,P − QED is given by [44,52] with the electron mass m e and the physical electromagnetic coupling constant e. Ψ and A µ are the fermion and the gauge boson fields, respectively. In the single-particle coordinate,Ĥ QED incorporates the transverse center-of-mass (CM) motion, which is not necessarily zero. The advantage of using the singleparticle coordinate is that each particle in the Fock space can be treated on equal footing, and this facilitates dealing with (anti-) symmetry among identical particles when higher Fock sectors are taken into account [44]. The disadvantage is that the CM motion of the system is mixed up with intrinsic motion.
In order to overcome this disadvantage, we add a constraint termĤ to the light-front QED Hamiltonian to factorize out the transverse CM motion from the intrinsic motion, where the CM motion is then governed by a harmonic oscillator Moreover, by subtracting the zero-point energy 2b 2 and multiplying by a Lagrange multiplier λ L , we are able to shift the excited states of CM motion to higher energy and retain low-lying states that contain only the simplest (gaussian) state of CM motion times the intrinsic motion. A detailed discussion of the CM factorization within the BLFQ framework can be found in Refs. [45,66]. The basis state used in BLFQ for a single particle is the direct product of the momentum eigenstate in the longitudinal direction, the two-dimensional harmonic oscillator basis state in the transverse direction and the helicity eigenstate. This choice is motivated by the success of ab initio calculation of nuclei structure [67,68] and also holographic QCD [46,69].
In the longitudinal direction, we impose (anti) periodic boundary condition on (fermions) bosons inside a box of length 2L, resulting in where k i is an integer (boson) or half-integer (fermion) 2 .
In the transverse direction, the basis states are the eigen- 2 We neglect the zero mode for bosons. For convenience we subtract 1 2 from state labels k i for fermions so that the longitudinal quantum number of all particles and their sum K will be designated as integers.
states of the two dimensional harmonic oscillator (2D-HO) Hamiltonian with two quantum numbers n i , m i and corresponding eigenenergies 2(2n i +|m i |+1)b 2 [66], where b is the HO scale parameter. Here n i represents the radial quantum number and m i the quantum number of orbital motion in the transverse plane. In the momentum space, the transverse basis states are given by Here, θ = arg(p ⊥ ), ρ = |p ⊥ |/b and L |m| n is the associated Laguerre function. Each single-particle basis state is identified using four quantum numbers: where λ i represents the light-front helicity [70] of the particle.
The basis states of a Fock sector with multiple particles are then expressed as the direct product of each singleparticle basis state in the Fock sector All the basis states have a well defined value of the total angular momentum projection This follows from the commutation betweenĴ Z and P − QED 3 . In addition, they all have the same total longitudinal momentum P + = i p + i . One then parameterizes P + using a dimensionless variable K = i k i such that P + = 2π L K. For a given particle i, the longitudinal momentum fraction is then represented as x i = p + i /P + = k i /K. Another element in reducing the infinite-dimensional basis to a finite-dimensional basis is the truncation of quantum numbers in the transverse direction. We retain the basis states satisfying: Note that N max not only restricts the HO quanta, but also introduces an IR (∼ b/ √ N max ) and UV (∼ b √ N max ) cutoff in the transverse direction. [52] Overall, there are three basis parameters in our calculation, the scale of the 2D-HO b, the transverse truncation N max , and the longitudinal resolution K. We will discuss the dependence of our results upon these parameters while discussing the numerical results in section IV.
A. ∆me and Z2 As mentioned before, we need to perform mass renormalization during the diagonalization and we employ a Fock sector dependent renormalization procedure [72,73]. The electron mass in the first Fock sector equals the bare electron mass, m e e = m 0 , and it depends on the regulators (cut-off) in the renormalization process. We numerically diagonalize the Hamiltonian matrix in an iterative scheme, where we adjust m e e in the first Fock sector | e to obtain the mass of the physical electron M e = 0.511 MeV as the mass of the entire system. Since in our Fock space truncation, Eq. (3), the photon cannot fluctuate into electron-positron pairs, for basis states in |eγ sector the electron mass remains the same as the physical value, m eγ e = M e . We introduce the mass counter-term ∆m e = m 0 − M e that represents the mass correction due to the quantum fluctuations to higher Fock sectors. According to the light-front perturbation calculation [64], ∆m e should increase with increasing basis truncation parameters, which leads to larger UV cutoff. This is confirmed in Fig. 1.
Furthermore, following the discussions in Refs. [43], [44], [50], [74], we rescale the naive TMDs by Z 2 to correct an artifact on the normalization of the LFWF in the | eγ Fock sector, which is introduced by our Fock sector truncation. Here, prime designates that the results are calculated directly by using the LFWF from Eq. (14). Z 2 is the electron wave function renormalization factor which in our truncation can be interpreted as the probability of finding a bare electron within a physical electron: The summation runs over all the basis states in the |e sector. Z 2 incorporates the contribution from the quantum fluctuation between the |e and |eγ sectors and goes to zero in the infinite basis limit. We show 1/Z 2 as a function of basis truncation parameters in Fig. 2. We comment that such sawtooth patterns are a familiar odd-even effect in the BLFQ formalism [44].

III. TMD IN QED
We treat the physical electron as a spin-1/2 composite particle, with a bare electron and a photon emerging from quantum fluctuation as its partons. The TMDs of a fermion inside the physical electron are defined through the fermion-fermion correlator function as 4 P, S |Ψ(0)U(0, z)ΓΨ(z)| P, S | z + =0 , where x = k + P + , and k ⊥ stands for the relative transverse momentum of the bare electron. Γ is a Dirac γmatrix. For the leading-twist TMDs, Γ = γ + , γ + γ 5 or iσ j+ γ 5 , corresponding to unpolarized, longitudinally polarized or transversely polarized fermion inside the physical electron respectively. The gauge link U(0, z) guarantees the gauge invariance of this non-local operator. | P, S ( P, S |) stands for the initial (final) state of the physical electron with momentum P and spin four-vector S.
We choose the light-cone gauge A + = 0. We work in a frame where the momentum and spin four-vector of the physical electron are represented as P = (P + , It is easy to verify that this covariant spin vector satisfies P · S = 0 and also S 2 = −1 if we choose the three-dimensional spin vector to be S = (S 1 , S 2 , S 3 ) = (sin θ cos ϕ, sin θ sin ϕ, cos θ) . (18) The spin vector state can then be transformed from the light-front helicity state [70] using the following formulae [7] | P, S This transformation is crucial because the eigenstates we obtain in BLFQ are the states associated with the lightfront helicities Λ = ±(1/2) rather than the spin vectors S. There exist eight leading-twist TMDs which are defined as Here, i, j run in the transverse plane, and anti-symmetric tensor 12 = − 21 = 1. However, if we set the gauge link to unity, two of them, the Boer-Mulders function h ⊥e 1 [75] and the Sivers function f ⊥e 1T [76] reduce to zero because of their T-odd property. Furthermore, we find that h ⊥e 1T mainly comes from the D-wave component and the D-wave component is zero in our current Fock sector truncation. Thus, h ⊥e 1T is also zero. This is consistent with the result of the perturbative light-front calculation [55].
Therefore, under the approximation that the gauge link is the identity operator and our current Fock sector truncation, there exist five non-zero TMDs at leadingtwist. We are now able to evaluate those TMDs using the LFWF in the | eγ Fock sector (ψ Λ λe,λγ (x e , k ⊥ e , x γ , k ⊥ γ ) in Eq. (2)) as [7,55] where k ⊥ R = k 1 + ik 2 and we use the abbreviation We omit the argument of LFWFs ({x i , k ⊥ i }) and TMDs (x, k ⊥ ) for the sake of conciseness.
The unpolarized TMD f e 1 describes the distribution of an unpolarized bare electron in an unpolarized physical electron, while the helicity TMD g e 1L gives the information of the momentum distribution of a longitudinally polarized bare electron when the physical electron is also longitudinally polarized. f e 1 and g e 1L are obtained by the sum of and the difference between the squares of S-and P-wave components under our Fock sector truncation, respectively. On the other hand, the transversity TMD h e 1 describes the correlation between the transversely polarized bare electron and the transversely polarized physical electron. Being chiral-odd in nature, h e 1 involves a helicity-flip of the bare electron from initial to final state, which is accompanied by a helicity flip of the physical electron in the same direction. This TMD receives contributions from the partial waves with L z = ±1. The TMD g e 1T is the distribution of a longitudinally polarized bare electron in a transversely polarized physical electron, whereas h ⊥e 1L provides the momentum distribution for the transversely polarized bare electron in a longitudinally polarized physical electron.
Meanwhile, as a spin-1/2 composite system, the three leading-twist PDFs of the physical electron can be retrieved by integrating f e 1 (x, k ⊥ ), g e 1L (x, k ⊥ ), and h e 1 (x, k ⊥ ) over k ⊥ ,

A. Averaging Method and Convergence Behaviour
In Fig. 3, we show the TMD f e 1 (x, k ⊥ ) in x and k ⊥ directions for fixed values of k ⊥ and x respectively. In these plots, the BLFQ computations are compared with those from the perturbative calculations in Ref. [55]. It can be noticed that in the transverse direction, the BLFQ results oscillate around the perturbative results, which can be viewed as a proxy for the experiment data in QED.
In this section we discuss the method we adopt to reduce the oscillation in the transverse direction i.e. with respect to k ⊥ and also study the convergence behaviour of our results as a function of the basis truncation parameters.
The oscillation in the TMD computed by the BLFQ method finds its root in the oscillating behaviour of the HO basis states that we use in the transverse plane. As mentioned in section II, N max , the truncation parameter in the transverse direction, together with b (the scale of the transverse 2D-HO basis), determines our transverse basis states. This motivates us to average over the BLFQ results at different N max and b to reduce these finite basis artifacts. We also have another truncation parameter K in the longitudinal direction. N max and K are truncation parameters introduced by our calculation method, so their effects must be studied to differentiate between truncation artifacts and the underlying physics extracted from the BLFQ results.
In our calculations, we set the 2D-HO scale b the same as the mass of physical electron i.e. b = M e . Here, we only discuss the dependence of our results on N max and K. In Fig. 3, we illustrate the BLFQ computations at different N max with K = 50 and b = M e . We see that if we only change N max , the phase of oscillation changes. The phase shift between computations at N max = n + 2 and at N max = n is nearly a half period. This property is very useful for approximately removing these truncation artifacts. In Fig. 4 we illustrate averages between the BLFQ computations at N max = n and N max = n + 2, with K = n. We further observe that the results obtained from averaging between {N max = n, n + 2} and between {N max = n + 2, n + 4} also shift from each other by nearly a phase of π. Based on these observations, we propose a two-step averaging method involving the BLFQ computations at {N max = n, n + 2, n + 4}, with K = n: 1. First take averages between the BLFQ computations at N max = n and N max = n + 2, with K = n. Then take averages between the BLFQ computations at N max = n + 2 and N max = n + 4, with also K = n.
2. The final results are obtained by averaging between those two averages produced in the first step.
The final averages for n = 50 are also shown in Fig. 4 and we see that they are very close to the perturbative results. In Fig. 5 we plot the BLFQ results after implementing the two-step averaging method using the same N max set but different K. We observe that the BLFQ results change minimally when K changes significantly since it is difficult to distinguish the difference between those curves at different K in Fig. 5. This suggests that our TMD results converge very quickly with longitudinal truncation parameter K. For simplicity and resolution in the longitudinal direction we from now on set K = n while implementing our averaging method.
In Fig. 6 we plot the BLFQ results averaged over N max = {n, n + 2, n + 4}, K = n for different n in the IR region (left panel) and also in the UV region (right panel) in the transverse direction. It is clear from the left panel of Fig. 6 that, after implementing our averaging method, the BLFQ results approach the perturbative result as N max increases. When we plot the BLFQ results in the UV region (right panel of Fig. 6), we observe that our result drops to zero abruptly after the UV cutoff (∼ b √ N max ). Since the UV cutoff increases as N max increases the agreement extends to the relevant higher momentum region.

B. TMDs
We evaluate the TMDs following Eqs. (23)- (27) and then obtain the BLFQ results by averaging over the BLFQ computations at N max = {100, 102, 104}, K = 100 and b = M e . If not stated otherwise, the BLFQ results that we show in the following are all obtained in the same way. In the left panel of Fig. 7, we show our five leading-twist TMDs as functions of k ⊥ for different values of x, and in the right panel of Fig. 7, we show them as functions of x for different values of k ⊥ . In these figures, we compare the BLFQ results with the perturbative results evaluated in Ref. [55]. We find that all the five non-zero TMDs obtained by the BLFQ approach are in excellent agreement with the perturbative results over the entire region of x and with moderately large |k ⊥ |. We argue that the discrepancy in the small |k ⊥ | region is due to relatively small truncation parameters leading to larger IR cutoffs, and the agreement improves when we relax the truncation parameters (see left panel of Fig. 6).    The three-dimensional structures of our BLFQ results after averaging are shown in Fig. 8. It is clearly seen there that all five non-zero TMDs have peaks near (x, k ⊥ ) = (1, 0). These peaks indicate the dominant probability is to find a bare electron in the physical electron in each of the different polarization configurations. For f e 1 , g e 1L and h e 1 , which appear very similar by eye, the peak in the transverse direction runs higher with increasing x but falls rapidly to zero or even negative region (for g e 1L ) after x decreases below 0.7. For g e 1T and h ⊥e 1L , the peak occurs very close to k ⊥ = 0 and falls rapidly towards zero as (k ⊥ ) 2 increases towards 0.01 (MeV 2 ).
PDFs of the bare electron within the physical electron are illustrated in Fig. 9, where again we find a good agreement between our BLFQ results after averaging and the perturbative results. We remark that in Fig. 9, those three PDFs are very close and in the limit x → 1 they basically become the same. This feature is also clearly seen if you take the same limit for the analytic perturbative equations, i.e. Eqs. (38), (39), (43) in Ref. [55]. In that limit, those three TMDs actually equal each other and thus their integration in the transverse direction also become the same.

C. Spin Densities in Momentum Space
By designating the three-dimensional spin vectors of the bare electron and physical electron by s = (s 1 , s 2 , s 3 ) and S = (S 1 , S 2 , S 3 ) respectively, we can define the bare electron's momentum distribution inside a physical electron with different polarization combinations via those five non-zero leading-twist TMDs as [55] ρ(x, k ⊥ , s, S) = Using this formula, we find that the momentum distribution for both the bare and the physical electron polarized in the same transverse direction can be computed by (f e 1 + h e 1 )/2. Other polarization combinations are (f e 1 + g e 1L )/2 for both polarized in the same longitudinal direction, (f e 1 + k y Me g e 1T )/2 for longitudinally polarized bare electron inside y−axis polarized physical electron and (f e 1 + k y Me h ⊥e 1L )/2 for the converse arrangement. In Figs. 10, 11, and 12, we present the bare electron densities in the transverse momentum plane at different values of x when both the bare and the physical electrons are unpolarized, longitudinally, and transversely polarized, respectively. We find that those three distributions in momentum space are azimuthally symmetric, having the maximum at the centre and rapidly fall off towards the periphery of the transverse momentum plane. Note that h ⊥e 1T and h ⊥e 1 vanish in our calculation which makes the transversely polarized density symmetric. In principle, with nontrivial gauge link they will not vanish and the transversely polarized density may not be symmetric. The densities grow larger and larger with increasing x, which is consistent with the right panels of Fig. 7. It is also observed from the right panels of Fig. 7 that f e 1 , g e 1L and h e 1 follow similar patterns in the longitudinal direction, though they are not identical. Thus, changes in the longitudinal direction in Figs. 10, 11, and 12 are also quite similar. However, longitudinally and transversely polarized distributions (Figs. 11, 12) have smaller values in the center with respect to unpolarized distribution (Fig. 10), owing to averaging f e 1 with g e 1L and f e 1 with h e  The distributions (f e 1 + k y Me g e 1T )/2 and (f e 1 + k y Me h ⊥e 1L )/2 are shown in Figs. 13 and 14, respectively, where we clearly see the distortion introduced by different polarization configurations of the bare and the physical electron. Both densities feature a significant dipole deformation arising from the terms k y Me g e 1T and k y Me h e⊥ 1L in y−axis. It can also be noticed that the distortion of a longitudinally polarized bare electron in a transversely polarized parent electron is opposite to that of a transversely polarized bare electron in a longitudinally polarized physical electron. The reason is that h ⊥e 1L is positive, while g e 1T is negative. Since both the TMDs have the peaks at large x, the distortions are stronger near x ∼ 0.8 compared to other values of x.
The oscillation of the BLFQ results after averaging (lower rows) around the perturbative results (upper rows) is also observed in Fig. 10 to Fig. 14. Although we utilize BLFQ results after averaging to calculate those density plots, we could still see the oscillation clearly. For Figs. 10, 11 and 12, the oscillation comes directly from the oscillation of f e 1 g e 1L and h e 1 in the small k ⊥ region even after averaging, which can be seen from the left panels of Fig. 7. And for Figs. 13 and 14, the original oscillation shown in Fig. 7 is small but are amplified by the transverse momentum k y in the definition of those distributions. We anticipate that BLFQ calculations in larger basis spaces, accompanied by a larger averaging domain, would further reduce these finite basis artifacts.

V. SUMMARY
In this paper, we have investigated the leading-twist TMDs of the physical electron from its light-front wave functions in the framework of Basis Light-front Quantization. These wave functions were obtained from the eigenvectors of the light-front QED Hamiltonian in the light-cone gauge. Our non-perturbative results have been compared with leading-order perturbative calculations. Our current study of the TMDs is performed within the lowest nontrivial Fock sectors. With a proper renormalization procedure and a rescaling of the naive TMDs correcting the artifacts introduced by the Fock space truncation, the BLFQ results are consistent with the perturbative calculations. During the calculation, we have also introduced an averaging method for reducing finite basis artifacts of our BLFQ results. This method may be used in the following non-perturbative investigation. In this study, the gauge link has been set to unity which leaves us five non-zero TMDs out of the eight leading-twist TMDs. Employing those non-zero TMDs, we have also studied the spin densities in the transverse momentum plane of a bare electron inside a physical electron with different polarization configurations. The quantitative results of the densities are also consistent with the perturbative results. These calculations constitute a comprehensive and accurate test of the BLFQ approach. For further investigation, future developments will focus on the inclusion of a nontrivial gauge link that will provide a prediction of the Boer-Mulders function h ⊥e 1 and the Sivers function f ⊥e 1T in QED. The main purpose of this study is to establish the foundation for studying the TMDs for strongly interacting systems in QCD, like the baryon and meson which are highly non-perturbative.