Particle distribution in intense fields in a light-front Hamiltonian approach

TL;DR

This paper employs the tBLFQ framework to study the non-perturbative real-time evolution of an electron in intense electromagnetic fields generated by ultra-relativistic nuclei, revealing significant higher-order effects for heavy nuclei.

Contribution

It demonstrates the application of the tBLFQ approach to realistic strong coupling QED problems involving intense fields from ultra-relativistic nuclei, highlighting non-perturbative effects.

Findings

Agreement with light-front perturbation theory for small nuclei

Identification of higher-order contributions for heavy nuclei

Observation of significant non-perturbative corrections in electron momentum evolution

Abstract

We study the real-time evolution of an electron influenced by intense electromagnetic fields using the time-dependent basis light-front quantization (tBLFQ) framework. We focus on demonstrating the non-perturbative feature of the tBLFQ approach through a realistic application of the strong coupling QED problem, in which the electromagnetic fields are generated by an ultra-relativistic nucleus. We calculate transitions of an electron influenced by such electromagnetic fields and we show agreement with light-front perturbation theory when the atomic number of the nucleus is small. We compare tBLFQ simulations with perturbative calculations for nuclei with different atomic numbers, and obtain the significant higher-order contributions for heavy nuclei. The simulated real-time evolution of the momentum distribution of an electron evolving inside the strong electromagnetic fields exhibits significant non-perturbative corrections comparing to light-front perturbation theory calculations. The formalism used in this investigation can be extended to QCD problems in heavy ion collisions and electron ion collisions.

Figures

Tables

Table 1. Table 1: Spinor background field potential vector contraction for different helicity configurations of the incoming electron (“1”) and the outgoing electron (“2”). We define the complex momentum as p ¯ = p x + i ​ p y ¯ 𝑝 superscript 𝑝 𝑥 𝑖 superscript 𝑝 𝑦 \bar{p}=p^{x}+ip^{y} .

(${\lambda }_2$,${\lambda }_1$)	$\overline{u}(p_2,{\lambda }_2){\gamma }^{\mu }u(p_1,{\lambda }_1){𝒜}_{\mu }\left(k\right)$
$\uparrow \uparrow$	$\frac{Z\left|e\right|}{ϵ\left(e^{-2y}{\left(k^{+}\right)}^2+k_{⟂}^2\right)}\left(2e^{-2y}+\frac{{\overline{p}}_1{\overline{k}}^{\ast }}{p_1^{+}{k}^{+}}+\frac{{\overline{p}}_1^{\ast }\overline{k}}{p_2^{+}{k}^{+}}\right)$
$\uparrow \downarrow$	$\frac{Z\left|e\right|}{ϵ\left(e^{-2y}{\left(k^{+}\right)}^2+k_{⟂}^2\right)}\frac{m_e{\overline{k}}^{\ast }}{k^{+}}\left(\frac{1}{p_1^{+}}-\frac{1}{p_2^{+}}\right)$
$\downarrow \uparrow$	$\frac{Z\left|e\right|}{ϵ\left(e^{-2y}{\left(k^{+}\right)}^2+k_{⟂}^2\right)}\frac{m_e\overline{k}}{k^{+}}\left(\frac{1}{p_2^{+}}-\frac{1}{p_1^{+}}\right)$
$\downarrow \downarrow$	$\frac{Z\left|e\right|}{ϵ\left(e^{-2y}{\left(k^{+}\right)}^2+k_{⟂}^2\right)}\left(2e^{-2y}+\frac{{\overline{p}}_1^{\ast }\overline{k}}{p_1^{+}{k}^{+}}+\frac{{\overline{p}}_1{\overline{k}}^{\ast }}{p_2^{+}{k}^{+}}\right)$
(${\lambda }_2$,${\lambda }_1$)	$\overline{u}(p_2,{\lambda }_2){\gamma }^i{\gamma }^{+}{\gamma }^{j}u(p_1,{\lambda }_1){𝒜}_i\left(k_2\right){𝒜}_j\left(k_1\right)$
$\uparrow \uparrow$	${\left(\frac{Z\left|e\right|}{ϵ}\right)}^2\frac{2{\overline{k}}_1^{\ast }{\overline{k}}_2}{\left(e^{-2y}{\left(k_1^{+}\right)}^2+k_1^{⟂}{}_{}{}^2\right)\left(e^{-2y}{\left(k_2^{+}\right)}^2+k_2^{⟂}{}_{}{}^2\right)}$
$\uparrow \downarrow$	0
$\downarrow \uparrow$	0
$\downarrow \downarrow$	${\left(\frac{Z\left|e\right|}{ϵ}\right)}^2\frac{2{\overline{k}}_1^{\ast }{\overline{k}}_2}{\left(e^{-2y}{\left(k_1^{+}\right)}^2+k_1^{⟂}{}_{}{}^2\right)\left(e^{-2y}{\left(k_2^{+}\right)}^2+k_2^{⟂}{}_{}{}^2\right)}$

Equations

$$|\beta ⟩=\sum _{\alpha }|\alpha ⟩⟨\alpha |\beta ⟩,$$

$$|e_{\text{phys}}⟩=a|e⟩+b|e\gamma ⟩+c|e\gamma \gamma ⟩+d|ee\bar{e}⟩+\dots ,$$

$$i\frac{\partial }{\partial x^{+}}|\psi ;x^{+}⟩=\frac{1}{2}{P}^{-}(x^{+})|\psi ;x^{+}⟩,$$

$$|\psi ;x^{+}{⟩}_I$$

$$|\psi ;0{⟩}_I=\sum _{\alpha }|\alpha ⟩c_{\alpha }(0),$$

$$|\psi ;x^{+}{⟩}_I:=\sum _{\alpha }{c}_{\alpha }(x^{+})|\alpha ⟩,$$

$$c(x^{+})$$

$$\mathcal{T}_{\scriptscriptstyle+}\exp\bigg{(}-\frac{i}{2}\int\limits_{0}^{x^{\scriptscriptstyle+}}\mathcal{M}\bigg{)}\rightarrow\mathcal{T}_{\scriptscriptstyle+}\prod_{n}\big{[}1-\tfrac{i}{2}\mathcal{M}(x^{\scriptscriptstyle+}_{n})\delta x^{\scriptscriptstyle+}\big{]}\;.$$

$$\mathcal{L}=-\frac{1}{4}{F}_{\mu \nu }{F}^{\mu \nu }+\bar{\Psi }(i{\gamma }^{\mu }{D}_{\mu }-m_e)\Psi ,$$

$$P^{-}=$$

$$P^{-}=P_{\text{QED}}^{-}+\mathcal{V}(x^{+}),$$

$$\rho (z,t,{\mathbf{x}}^{\perp })=Z|e|\delta (z-vt)\delta ({\mathbf{x}}^{\perp }),$$

$$\mathbf{j}(z,t,{\mathbf{x}}^{\perp })=Z|e|v\hat{\mathbf{z}}\delta (z-vt)\delta ({\mathbf{x}}^{\perp }).$$

$$({\nabla }^2-{\partial }_t^2){\mathcal{A}}^0=-\rho ,$$

$$({\nabla }^2-{\partial }_t^2)\mathbf{A}=-\mathbf{j},$$

$${\mathcal{A}}^{-}=$$

$${\mathcal{A}}^i=$$

$$|{\varphi }_0⟩=|p_0^{+},G({\sigma }_0,{\mathbf{p}}_0^{\perp }),{\lambda }_0⟩.$$

$${\varphi }_0^{\perp }({\mathbf{p}}^{\perp })=⟨{\mathbf{p}}^{\perp }|G({\sigma }_0,{\mathbf{p}}_0^{\perp })⟩=\frac{1}{\sqrt{\pi }{\sigma }_0}{e}^{-\frac{({\mathbf{p}}^{\perp }-{\mathbf{p}}_0^{\perp }{)}^2}{2{\sigma }_0^2}}.$$

$$|{\varphi }_f⟩=|p_f^{+},{\mathbf{p}}_f^{\perp },{\lambda }_f⟩,$$

$$\int {\mathrm{d}}^2{\mathbf{p}}^{\perp }⟨p_f^{+},{\mathbf{p}}_f^{\perp },{\lambda }_f|{\mathcal{W}}_1|p_0^{+},{\mathbf{p}}^{\perp },{\lambda }_0⟩{\varphi }_0^{\perp }({\mathbf{p}}^{\perp }).$$

$$\sum _{\alpha ,{\alpha }^{\prime }}⟨p_f^{+},{\mathbf{p}}_f^{\perp },{\lambda }_f|{\alpha }^{\prime }⟩⟨{\alpha }^{\prime }|{\mathcal{W}}_1|\alpha ⟩⟨\alpha |p_0^{+},G({\sigma }_0,{\mathbf{p}}_0^{\perp }),{\lambda }_0⟩,$$

$$\sum _{n,m}{\tilde{\Phi }}_{nm}^{\ast }({\mathbf{p}}^{\perp }){\tilde{\Phi }}_{nm}({\mathbf{q}}^{\perp })=(2\pi {)}^2{\delta }^{(2)}({\mathbf{p}}^{\perp }-{\mathbf{q}}^{\perp }),$$

$$\int \frac{{\mathrm{d}}^2{\mathbf{p}}^{\perp }}{(2\pi {)}^2}\sum _{n,m}{\tilde{\Phi }}_{nm}^{\ast }({\mathbf{p}}^{\perp }){\tilde{\Phi }}_{nm}(0)$$

$$\mathcal{S}=1-\frac{i}{2}\underset{T\to \infty }{\lim }{\int }_0^{T}d{x}^{+}{P}_{int}^{-}(x^{+})+...$$

$${\Gamma }_T(i\to f)=\frac{P_T(i\to f)}{T},$$

$$\mathcal{T}_{\scriptscriptstyle+}e^{-\frac{i}{2}\int\limits_{0}^{x^{+}}P^{\scriptscriptstyle-}}\rightarrow\mathcal{T}_{\scriptscriptstyle+}\prod_{n}\big{[}1-\tfrac{i}{2}P^{\scriptscriptstyle-}(x^{+}_{n})\delta x^{+}\big{]}\;.$$

$$f_0(p^{\perp })=\frac{p^{\perp }}{2\pi }|{\tilde{\Phi }}_{00}^b({\mathbf{p}}^{\perp }){|}^2,$$

$$x^{\pm }=x^0\pm x^3,{\mathbf{x}}^{\perp }=(x^1,x^2)$$

$$H_{HO}^{2d}=\frac{p_{\perp }^2}{2M}+\frac{1}{2}M{\Omega }^2{x}_{\perp }^2,$$

$${\Phi }_{nm}^b(\rho ,\varphi )=(-1{)}^n{i}^{|m|}{f}_{nm}(\rho ){\chi }_m(\varphi ),$$

$$f_{nm}^b(\rho )=b\sqrt{2}\sqrt{\frac{n!}{(n+|m|)!}}{e}^{-b^2{\rho }^2/2}(b\rho {)}^{|m|}{L}_n^{|m|}(b^2{\rho }^2),$$

Full text

I Introduction

Real-time quantum field theory (QFT) in the non-perturbative regime is pivotal for understanding recent experimental discoveries in modern high energy nuclear colliding facilities, e.g., the Large Hadron Collider (LHC) and the Relativistic Heavy Ion Collider (RHIC). For instance, thorough understanding of jet quenching and heavy quarkonium suppression in heavy ion collisions requires detailed knowledge of how quarks and gluons interact with the evolving hot medium created by the colliding nuclei Qin:2015srf ; Rapp:2008qc . Other examples of time-dependent non-perturbative problems are, QED in ultra-intense laser fields Heinzl:2008an ; DiPiazza:2011tq and QCD in strong magnetic fields Skokov:2009qp ; Tuchin:2013ie .

Stationary state QFT problems in the non-perturbative regime are themselves challenging, and basis light-front quantization (BLFQ) has emerged as a promising framework to solve non-perturbative QFT eigenstates from first-principles Honkanen:2010rc ; Vary:2009gt ; Vary:2016emi . By employing the light-front Hamiltonian formalism, BLFQ enjoys advantages of light-front dynamics and of non-perturbative quantum many-body theory at the same time. It shares many advantageous features with discretized light-front quantization (DLCQ) Pauli:1985ps and ab initio nuclear structure calculations, e.g., the no-core shell model (NCSM) Navratil:2000ww ; Navratil:2000gs ; Barrett:2013nh . Additionally, the basis approach explicitly retains kinematical symmetries of the system and could lead to significant reduction of numerical efforts. In addition to providing the mass eigenstates, the light-front Hamiltonian formalism generates the associated light-front amplitudes which can then be applied to determine spin structures, electromagnetic form factors and generalized parton distributions of hadrons and other observables. BLFQ has been successfully applied to QED problems, e.g., the electron anomalous magnetic moment Honkanen:2010rc ; Zhao:2014xaa , the positronium system Wiecki:2014ola and GPDs of the electron and strong coupling positronium Chakrabarti:2014cwa ; Adhikari:2016idg . Recently, BLFQ was applied to the heavy quarkonium system with a confinement potential inspired from anti-de Sitter/conformal field theory (AdS/CFT) along with the one-gluon exchange interaction from the QCD Hamiltonian Li:2015zda . The spectroscopy and decay constants obtained from the BLFQ approach are comparable to experimental measurements and other established methods.

Distinct from the Lagrangian formulations, the Hamiltonian approach permits access to the real-time evolution of quantum states. Therefore, time-dependent basis light-front quantization (tBLFQ) is a natural extension of the BLFQ formalism. The tBLFQ formalism has been successfully applied to nonlinear Compton scattering by employing a simple anzatz for the time-dependent and intense laser field as a classical background field. Zhao et al.  have illustrated that tBLFQ enables real-time accessibility to intermediate quantum states of the electron by showing the evolution of the invariant mass of the electron and photon Fock state, which agrees with a perturbative calculation in the small coupling limit Zhao:2013jia ; Zhao:2013cma .

It is well known that a classical description can capture the most substantial physics of the system when the occupation number is large in quantum phase space. One example is laser physics Sarachik:1970ap . Another renowned example is the Color Glass Condensate (CGC) McLerran:1993ni ; JalilianMarian:1997jx ; JalilianMarian:1997gr ; Iancu:2000hn , a classical effective theory of QCD, where the small-$x$ partons are treated as classical fields generated by large-$x$ partons. Quantum effects are then treated as higher-order corrections to the classical calculations. The universal description of saturated gluons in hadrons based on CGC effective theory is able to explain a wide range of phenomena in deep-inelastic scattering and hadron-hadron collision experiments at high energies. For a recent review on this topic, see Ref. Gelis:2010nm . We foresee tBLFQ to be a very useful tool to study high-energy heavy ion collisions and electron ion collisions in conjunction with CGC effective theory.

As a first realistic application of the tBLFQ framework to high-energy heavy ion collisions, we investigate QED effects in which the role of strong electromagnetic (e.m.) field is yet to be understood quantitatively Tuchin:2012mf ; Tuchin:2013apa ; Tuchin:2013ie . To be more specific, we study features of realistic electromagnetic fields generated by an ultra-relativistic heavy ion using the tBLFQ formalism. We solve the time-evolution of the quantized field of an electron inside such classical fields. The coupling between the electron and the strong e.m. field is at order $Z\alpha_{\text{em}}$ with $Z$ the atomic number of the nucleus and $\alpha_{\text{em}}$ the electromagnetic coupling constant. A non-perturbative approach is preferred for the strong coupling QED problem when $Z$ is large. For instance, the coupling between electron and e.m. field generated by gold nucleus is around $0.6$. We focus on demonstrating the non-perturbative features of the tBLFQ framework. This investigation also serves as a stepping-stone for future applications of tBLFQ to QCD problems in high energy nuclear collisions.

The paper is organized as follows. In section II, we provide a brief review on the background of this investigation. Then in section III we compare the tBLFQ simulation to light-front perturbation theory (LFPT). Effects on physical observables of high energy nuclear experiments are shown in section IV. In section V, we summarize our results and discuss additional applications of the tBLFQ framework in heavy ion collisions and electron ion collisions.

II Background

First we briefly review the key features of BLFQ and tBLFQ, and the components of the QED Hamiltonian with classical background fields. We refer readers to Ref. Zhao:2013cma for details. We will close this section by discussing some general properties of the intense electromagnetic fields generated by an ultra-relativistic heavy ion.

II.1 Basis Light-front Quantization

Obtaining the invariant mass eigenstates in a light-front Hamiltonian matrix approach has shown significant promise Brodsky:1997de ; Wiecki:2014ola ; Li:2015zda ; Vary:2016emi . The primary advantage of BLFQ is that, by adopting a basis with the same symmetries of the system under investigation, we can reduce the numerical efforts required for an accurate representation of the Hamiltonian.

The choice of basis is arbitrary as long as it is orthogonal and complete. One of the many choices is the two dimensional harmonic oscillator (‘2D-HO’) basis in the transverse direction and the discretized plane-wave basis in the longitudinal direction. Each single-particle basis state can be identified using four quantum numbers, $\bar{\alpha}=\{k,n,m,\lambda\}$. The longitudinal momentum of the particle is characterized by the first quantum number $k$. In the longitudinal direction $x^{\scriptscriptstyle-}$, we constrain the system to a box of length $2L$, and impose (anti-) periodic boundary conditions on (fermions) bosons. As a result, the longitudinal momentum $p^{\scriptscriptstyle+}=2\pi k/L$ is discretized, where the dimensionless quantity $k=1,2,3,...$ for bosons and $k=\frac{1}{2},\frac{3}{2},\frac{5}{2},...$ for fermions. We have neglected zero modes for bosons. The length parameter $L$ should be chosen to cover the longitudinal extent of the system, we will discuss it in section II.4. The next two quantum numbers, $n$ and $m$, depict radial excitation and angular momentum projection, respectively, of the particle within the 2D-HO basis in the transverse direction. The 2D-HO basis may be defined by two parameters, mass $M$ and frequency $\Omega$. However, we adopt a single HO parameter $b:=\sqrt{M\Omega}$, since our transverse modes depend only on $b$ rather than on $M$ and $\Omega$ individually. The state carrying quantum number $n$ and $m$ has HO eigenenergy $E_{n,m}=(2n+|m|+1)\Omega$, see Appendix A for details.

The many-particle basis states $|\,\alpha\,\rangle$ in each Fock sector are direct products of single-particle states. Such basis was initially designed for the QCD bound state problem and was supported by AdS/CFT correspondence with QCD. It has been shown, such a choice of basis allows one to encode the following three symmetries of QED Hamiltonian. First, translational symmetry in the $x^{\scriptscriptstyle-}$ direction, i.e. conservation of longitudinal momentum $P^{\scriptscriptstyle+}$. Second, rotational symmetry in the transverse plane i.e. the longituidinal projection of angular momentum $J^{3}=J^{3}_{\text{orbital}}+J^{3}_{\text{spin}}$ is conserved. Finally, lepton number conservation, i.e. net fermion number is conserved and so is the total charge. Therefore the eigenspace of QED can be grouped into segments with definite eigenvalues ${K,M_{j},N_{f}}$ for the operators ${P^{\scriptscriptstyle+},J^{3},Q}$, respectively.

The physical QED eigenstates, written as $|\,\beta\,\rangle$, are represented as the superposition of the basis states,

[TABLE]

with eigenstates and basis states belonging to the same segment. Coefficients $\langle\,\beta\,|\,\alpha\rangle$ are obtained by diagonalizing $P^{\scriptscriptstyle-}_{\text{QED}}$ in the basis representation. To this end, we require a finite dimension representation of the QED Hamiltonian that is achieved through the following basis reduction procedures.

First, by taking into account the conserved quantities and selection rules, one determines which subset of basis states will contribute to a desired observable. For this observable, one needs to work in a finite number of segments of QED eigenspace without any information loss. Second, because of the fact that even a single segment has an infinite number of degrees of freedom, we need to truncate the basis which inevitably introduces loss of precision in our calculated observables. We implement two levels of truncation scheme as follows.

i) Fock-sector truncation. Take the physical electron state as an example. Schematically, a physical state with $N_{f}=1$ has the following Fock sector expansion

[TABLE]

containing the bare electron $|\,e\,\rangle$ and its dressed states $|\,e\gamma\,\rangle$, $|\,e\gamma\gamma\,\rangle$ and $|ee\bar{e}\rangle$ etc. We explicitly assume that higher Fock-sectors give insignificant contributions to the low-lying eigenstates in which we are mostly interested, with an appropriate renormalization procedure implemented. Such assumption is motivated by the success of perturbation calculations in QFT. Furthermore, the dominance of contributions to physical observables from lower Fock sectors has been shown in scalar Yukawa model even in the non-perturbative regime Li:2015iaw , which indirectly justifies our Fock sector truncation scheme in QED. In the following calculations we make the simplest truncation, by keeping only the single electron Fock sector in Eq. (2), since the interactions between the fermion and the photon are suppressed by $1/Z$ comparing to the interaction between the fermion and the classical field generated by the nucleus, where $Z$ is the atomic number of the nucleus. We leave the corrections from higher Fock sectors to a future study.

ii) Truncation within Fock-sectors. Within each Fock-sector, further truncations are still needed to reduce the basis to a finite dimension. As mentioned, we impose (anti-) periodic boundary conditions on (fermions) bosons in a longitudinal box with length $2L$. Consequently, the longitudinal momentum $p^{\scriptscriptstyle+}$ of single particles can only take discrete values. We then introduce a truncation parameter $K$ on the longitudinal direction such that, $\sum_{l}k_{l}\leq K$, where $k_{l}$ is the longitudinal momentum quantum number of $l$-th particle in the basis state. Note that systems with larger $K$ have simultaneously higher ultra-violet (UV) and lower infra-red (IR) cutoffs in the longitudinal direction. In the transverse direction, we require the total transverse quantum number $N_{\alpha}=\sum_{l}(2n_{l}+|m_{l}|+1)$ for multiparticle basis state $|\,\alpha\,\rangle$ satisfies $N_{\alpha}\leq N_{\text{max}}$, where $N_{\text{max}}$ is a chosen truncation parameter.

We thus attain a finite dimensional representation of the QED Hamiltonian in the BLFQ basis. The continuum limit can be achieved by extending $K$ and $N_{\text{max}}$ to infinity. The dependence on the parameter $L$ should be weak as long as it covers the longitudinal extent of the system.

II.2 Time-dependent Basis Light-front Quantization

The state of a quantum system at a later time is related to its state at an earlier time by the Schrödinger equation, which takes the following form,

[TABLE]

in light-front dynamics. The Schrödinger equation can be solved in either the interaction picture or the Schrödinger picture. Physical observables should not depend on the pictures we employed for the time evolution. However, for a particular problem, one picture may be numerically advantageous over another. For instance, if the Hamiltonian has a non-trivial time dependence, working in the Schrödinger picture may be more numerically efficient since the interaction picture would require calculating the Hamiltonian in the physical eigenstates at every time step. On the other hand, if the interaction is much smaller than the kinetic Hamiltonian, then a finer time step is required in the Schrödinger picture to produce the same precision as in the interaction picture. Of course the choice also depends on the physical observable of interest, e.g., the interaction picture could be more capable in describing bound states. In a word, the time evolution picture should be chosen according to the problem itself. In this investigation, since we are interested in the effects of external fields, we work in the interaction picture for time evolution. Its formal solution is,

[TABLE]

where $\mathcal{T}_{\scriptscriptstyle+}$ is light-front time ordering operator and $V_{I}$ is the interaction Hamiltonian in the interaction picture, with the subscript $I$ indicates the interaction picture. We can expand the initial state in the BLFQ basis,

[TABLE]

where $c_{\alpha}(0)\equiv\langle\,\alpha\,|\,\psi;0\rangle_{I}$. The coefficients of the state at later times can be expanded as

[TABLE]

in the BLFQ basis. Its coefficients will be solved through,

[TABLE]

where $\mathcal{M}$ is a finite dimensional representation of the Hamiltonian operator in BLFQ basis, $\mathcal{M}_{\alpha\alpha^{\prime}}=\langle\,\alpha\,|V_{I}|\,\alpha^{\prime}\,\rangle$. The time-evolution operator then is decomposed into small steps in light-front time $x^{+}$, with step size $\delta x^{+}$,

[TABLE]

A higher-order difference scheme Askar:1978 ; Iitaka:1994 is implemented to ensure numerical stability and precision, refer to Appendix B for details. The continuum limit corresponds to the limit taking step size $\delta x^{\scriptscriptstyle+}\rightarrow 0$.

II.3 The Light-Front QED Hamiltonian

Starting from the QED Lagrangian with an additional background field,

[TABLE]

where $D_{\mu}\equiv\partial_{\mu}+ieC_{\mu}$ and $C_{\mu}=\mathcal{A}_{\mu}+A_{\mu}$ is the sum of the background and quantum gauge fields, respectively. In this paper, $\mathcal{A}_{\mu}$ is the electromagnetic field generated by the nucleus with atomic number $Z$. Note that $F_{\mu\nu}$ is calculated from $A_{\mu}$ alone. Working in the light-front gauge, the full Hamiltonian is then derived as Brodsky:1997de ; Zhao:2013cma ,

[TABLE]

The first three lines are the QED light-front Hamiltonian, $P^{-}_{\text{QED}}$. In order, each of the first five terms in Eq. (II.3) represents the fermion kinetic energy $T_{f}$, photon kinetic energy $T_{\gamma}$, vertex interaction $W_{1}$, instantaneous-photon interaction $W_{2}$ and instantaneous-fermion interaction $W_{3}$ respectively. The last two lines contain the four new interactions generated by the classical background field $\mathcal{A}$, we label them as $\mathcal{W}_{1}$, $\mathcal{W}_{2}$, $\mathcal{W}_{3}$ and $\mathcal{W}_{4}$ respectively. Since we only keep the leading Fock-sector which contains one single fermion, only $T_{f}$, $\mathcal{W}_{1}$ and $\mathcal{W}_{2}$ enter our calculation. In the following, we only keep the relevant terms in the QED Hamiltonian,

[TABLE]

where $P^{-}_{\text{QED}}=T_{f}$ and $\mathcal{V}(x^{\scriptscriptstyle+})=\mathcal{W}_{1}(x^{\scriptscriptstyle+})+\mathcal{W}_{2}(x^{\scriptscriptstyle+})$ throughout our discussion here. Note $\mathcal{W}_{1}$ is first order in $Z\alpha_{\text{em}}$ and $\mathcal{W}_{2}$ is second order in $Z\alpha_{\text{em}}$, with $\alpha_{\text{em}}\sim 1/137$ being the electromagnetic coupling constant.

Our particular truncation of Fock-sector also simplifies the problem as we can take physical values for the electron mass and charge. If one works with higher Fock-sectors with both electron and photon, proper renormalization is required. One feasible choice is the sector dependent scheme Karmanov:2008br ; Karmanov:2012aj ; Chabysheva:2009ez , which has been applied to the QED Hamiltonian when calculating the electron anomalous magnetic moment, for which the result agrees with the Schwinger value within 1% Zhao:2014xaa .

II.4 Electromagnetic Fields Generated by Relativistic Heavy Ion

The charge densities and current densities of one ion with atomic number $Z$ moving along the $z$-axis with velocity $v$ are

[TABLE]

The four vector potential of the fields obeys,

[TABLE]

where we omit the vacuum permittivity and permeability in natural units. In the light-cone gauge, in terms of $k^{\scriptscriptstyle+},{\bm{k}}^{\perp},x^{\scriptscriptstyle+}$, the solutions of the above equations are,

[TABLE]

where $y=\frac{1}{2}\ln(P^{\scriptscriptstyle+}/P^{\scriptscriptstyle-})$ is the rapidity of the heavy ion, with $P^{\mu}$ the momentum four-vector of the heavy ion.

Let us discuss the spatial distribution of the potential before we proceed with our calculation. In modern high energy collision facilities, particles are accelerated to the ultra-relativistic regime. For example, at RHIC, the center of mass energy of the collisions reaches $200$ GeV per nucleon, the rapidity of the colliding nuclei is $y\approx 5.3$. The energy is even higher at the LHC, e.g.  the rapidity of the colliding particles are $y\sim 9.5$ at center of mass energy $1.4$ TeV. Thus $e^{-2y}$ is a small (large) number when the heavy ion is moving along positive (negative) $z$-axis. Apparently, the potential generated by a heavy ion moving along the positive $z$-axis is almost stationary with a period $2\pi e^{2|y|}/k^{\scriptscriptstyle+}$ in $x^{\scriptscriptstyle+}$, and it has a very narrow extent in the longitudinal direction, see Fig. 1 . Contrarily, the potential generated by a heavy ion moving along negative $z$-axis is oscillating rapidly with a period $2\pi e^{-2|y|}/k^{\scriptscriptstyle+}$ in $x^{\scriptscriptstyle+}$, while it has a very wide extent in the longitudinal direction, see Fig. 1.

Although the potential profiles for a heavy ion moving along the positive and negative direction are very different from a first look, physical observables for the same process must be independent of such mathematical treatment. In tBLFQ, to achieve an accurate description of the same process, a larger coverage in $x^{\scriptscriptstyle+}$ ($x^{\scriptscriptstyle-}$) is necessary for a heavy ion moving in the positive (negative) $z$-direction and their continuum limits should be equivalent. In the following discussion, we assume the heavy ion is moving along the positive $z$-axis. First, it is easier to handle the time evolution of a quasi stationary potential. Second, the potential is concentrated in a smaller region in $x^{\scriptscriptstyle-}$. Consequently a moderate truncation parameter $L$ is sufficient to enclose the potential and cover a wide region in longitudinal momentum at the same time.

We are now ready to discuss the truncation parameter $L$ introduced in section II.1. The guidelines are that, the potential outside the box only makes inconsequential corrections to the process we are interested in, while at the same time a smaller $L$ provides larger longitudinal momentum coverage for the same truncation parameter $K$. Eqs. (14) clearly suggest that the longitudinal extent of the potential in $x^{\scriptscriptstyle-}$ depends on the transverse momentum we are interested in. We will specify our choice of $L$ for each calculation we perform in the following.

III Comparison to Light-Front Perturbation Theory

The equivalence of the LFPT and the covariant perturbation theory has been established decades ago Chang:1973qi . Furthermore, QED in the perturbative regime has been verified up to very high precision by various experiments Mohr:2008fa . Thus we can check the validity of our formalism using LFPT in the perturbative regime and study the numerical error introduced by the truncations and time-step discretization. In addition, we can study higher-order contributions by comparing tBLFQ simulations to LFPT calculations.

III.1 Comparison to Momentum Basis

The Hamiltonian matrix elements $\langle\,\alpha^{\prime}\,|\mathcal{V}|\,\alpha\,\rangle$ for the potential in Eq. (14) in the BLFQ basis can be calculated algebraically and the detailed expressions are presented in Appendix C. The first check would be that the transition amplitudes induced by the interaction $\mathcal{V}$ between particular initial and final states are consistent in the BLFQ basis and momentum basis, with sufficiently large $N_{\text{max}}$. We compare the leading order in the coupling between electron and background field $\alpha\equiv Z\alpha_{\text{em}}$ so that only $\mathcal{W}_{1}$ is relevant. In the longitudinal direction we adopt the discretized momentum basis, with which we approach the continuum longitudinal momentum limit when $K$ increases. For simplicity in the HO basis, we adopt a wave packet which is a Gaussian in the transverse direction. The Gaussian packet is centered at ${\bm{p}}_{0}^{\perp}$ and the width of the Gaussian wave packet $\sigma_{0}$ can be chosen independent of the 2D-HO parameter in the BLFQ basis. When $\sigma_{0}\rightarrow 0$ the initial state becomes a representation of the transverse momentum eigenstate. In the longitudinal direction we use its momentum eigenstate. Thus, the initial state is labeled by the following quantum numbers,

[TABLE]

Its normalized wavefunction in the transverse direction is,

[TABLE]

The transition amplitude $\langle\,\phi_{f}\,|\mathcal{W}_{1}|\,\phi_{i}\,\rangle$ to a momentum eigenstate,

[TABLE]

can be calculated by integrating over the initial transverse momentum distribution,

[TABLE]

It can also be calculated in the BLFQ basis as follows,

[TABLE]

where $\langle\,\alpha\,|\,p^{\scriptscriptstyle+}_{0},G(\sigma_{0},{\bm{p}}_{0}^{\perp}),\lambda_{0}\rangle$ can be calculated analytically. In principle, Eqs. (18) and (19) are identical if we sum over all $\alpha$ states. In practice, we can only perform calculations using the truncation scheme explained in section II.1. It is then necessary to check the behavior of the numerical uncertainty as we increase our truncation parameters.

We study the transition amplitude of an electron influenced by the field generated by a gold nucleus moving along the positive $z$-axis with rapidity $y=5.3$. We study transition amplitudes calculated in the BLFQ basis and the momentum basis from various initial states and final states, the comparisons show similar convergence behavior. As an example, we consider the initial state $|\,\phi_{0}\,\rangle$ at $x^{\scriptscriptstyle+}=0$ with $p^{\scriptscriptstyle+}_{0}=3\pi/10$ GeV, $\lambda_{0}=1/2$, and a Gaussian wave packet centered at ${\bm{p}}_{0}^{\perp}=(b/2,0)$ with width $\sigma_{0}=b/2$. We set the HO parameter in the BLFQ basis to $b=1000m_{e}=0.511$ GeV. We calculate the amplitude to the final state $|\,\phi_{f}\,\rangle$ which has the following quantum number, $p^{\scriptscriptstyle+}_{f}=7\pi/10$ GeV,$\lambda_{f}=1/2,{\bm{p}}_{f}^{\perp}=(p^{x}_{f},0)$ in the BLFQ basis with truncation $N_{\text{max}}$ and compare it to a momentum basis calculation. Fig. 2 shows the comparison between the BLFQ basis results with different $N_{\text{max}}$ and the momentum basis calculation. The difference between those two bases decreases as we increase $N_{\text{max}}$. However, even for $N_{\text{max}}=50$, we find their discrepancy is not negligible. On the other hand, if we average results of two BLFQ basis calculations with consecutive even $N_{\text{max}}$, the results show excellent agreement with the momentum basis calculation at $N_{\text{max}}=50$. This indicates averaging over $N_{\text{max}}$ can effectively enhance the convergence to momentum basis calculation.

Mathematically, the challenge in comparing HO basis and momentum basis calculation is rooted in the fact that HO functions are square-integrable while plane waves are not. Intuitively, the averaging procedure we adopted above can be illuminated by the following example. It is well known that a 2D delta function has its 2D-HO wavefunction representation,

[TABLE]

the equality is exact if $n$ and $m$ are summed over all possible values. However if we constrain the domain of $n$ and $m$ by requiring $N_{\text{total}}=2n+|m|+1\leq N_{\text{max}}$, the following integral

[TABLE]

is oscillating between the discrete results [math] and $2$ as $N_{\text{max}}$ increases. Such a behavior provides an argument for averaging over consecutive $N_{\text{max}}$ results when comparing HO basis calculations to momentum basis calculations.

We have shown that when $N_{\text{max}}$ and $K$ are sufficiently large, the BLFQ basis results exhibit good agreement with the momentum basis results. At this stage, we can now concentrate on the behavior of discretized numerical time evolution scheme where all calculations are performed in the BLFQ basis.

III.2 Transition Rate

Perturbation theory provides contributions from the interaction through a power series in the coupling constant $\alpha$, which is expected to converge for sufficiently small $\alpha$. Applying perturbation theory to the S-matrix, we have,

[TABLE]

where we have only kept the leading order (LO) contribution in $\alpha$. Contrary to the infinite time limit in perturbation theory, we can only evolve for a finite time using the discretized numerical time evolution scheme described in Eq. (7). Extrapolation to the infinite time limit requires infinite energy resolution which can be achieved only in the limit with $N_{\text{max}}$ and $K$ approaching infinity. With finite $N_{\text{max}}$ and $K$ we make a compromise by comparing to the perturbative calculation without taking the infinite time limit.

The most important physical observable measured in scattering experiments is the cross section, which is related to the transition rate by a flux factor. The transition rate is defined as the transition probability from state $|\,i\,\rangle$ to $|\,f\,\rangle$ divided by the time $T$ during which the interaction is active,

[TABLE]

where $P_{T}(i\rightarrow f)=|\langle\,f\,|\mathcal{S}|\,i\,\rangle|^{2}$. The scattering matrix can be calculated using either LFPT or the tBLFQ formalism. One advantage of the tBLFQ approach is that it does not rely on any expansion in the coupling constant of the interaction thus it is a legitimate approach for both weak and strong interactions. The applicability of tBLFQ for time-dependent non-perturbative problems is rooted in the numerical schemes we adopted for basis construction and time evolution. We have discussed extensively the non-perturbative feature of the BLFQ approach in section II.1. Therefore, we now focus on how non-perturbative effects are incorporated by our discretized time evolution scheme.

Conceptually, all numerical time evolution schemes are implemented by decomposing the time evolution operator into many small steps with step size ($\delta x^{\scriptscriptstyle+}$) in light-front time $x^{\scriptscriptstyle+}$,

[TABLE]

Contributions from up to $n$-th order in $\alpha$ are preserved, with $n$ the total number of time steps. The resummation up to all orders in $\alpha$ is achieved by taking the limit $\delta x^{\scriptscriptstyle+}\rightarrow 0$.

As a demonstration of the non-perturbative feature of the tBLFQ approach, we consider the transition of an electron between two physical QED eigenstates $|\,\beta\,\rangle$, as defined in Eq. (1), in the fields generated by different nuclei. We construct the BLFQ basis states as follows. We require the electron to be in the segments with definite eigenvalues $m=0$ and $\lambda=1/2$. We can make such a choice for the following two reasons. First, the potential is azimuthally symmetric, thus the total angular momentum projection $M_{j}=m+\lambda$ of the evolving electron is conserved. Second, the spin projection of the electron $\lambda$ is approximately conserved, since the helicity flip processes are suppressed by $m_{f}/p^{\perp}$ comparing to the helicity non-flip processes, where $m_{f}$ is the mass of the fermion and $p^{\perp}$ is its typical transverse momentum, see Table. 1 for detail. We set the longitudinal box length to be $L=10$ GeV${}^{-1}\approx 2$ fm. As discussed above, this is sufficient for the potential in our application. We set the HO parameter to be $b=1000m_{e}=0.511$ GeV, which is chosen as a representative of the typical transverse momentum of particles observed in heavy ion collisions. The BLFQ basis is constructed using $N_{\text{max}}=32$, $K=32$. The physical QED eigenstates are then obtained by diagonalizing the fermion kinetic energy in this BLFQ basis.

We take each nucleus to be moving along the positive $z$-axis with $y=5.3$. Since the nucleus is moving almost along $x^{\scriptscriptstyle+}$, we approximate the generated potential as a static potential during the time interval our calculation is performed. We then assume an electron enters such fields at time $x^{\scriptscriptstyle+}=0$ in a physical QED $|\,\beta_{i}\,\rangle$ of $P^{\scriptscriptstyle-}_{\text{QED}}$, with energy $P^{\scriptscriptstyle-}_{\beta,i}=0.455$ GeV, which belongs to the segment of BLFQ basis states with $k=43/2$. We calculate the transition rate of the electron to the final state $|\,\beta_{f}\,\rangle$ with energy $P^{\scriptscriptstyle-}_{\beta,f}=0.955$ GeV, which belongs to the segment of BLFQ basis states with $k=45/2$. We have only considered the transition induced by the $\mathcal{W}_{1}$ term of the Hamiltonian. In Fig. 3, we show the transition rate from the initial eigenstate $|\,\beta_{i}\,\rangle$ with energy $P^{\scriptscriptstyle-}_{\beta,i}=0.455$ GeV to the final state $|\,\beta_{f}\,\rangle$ with energy $P^{\scriptscriptstyle-}_{\beta,f}=0.955$ GeV, for nuclei with atomic number $Z=1,28,56,79$ as a function of exposure time T.

We make the following observations. First, in the small coupling regime, the tBLFQ calculation agrees with the NLO LFPT calculation, and both of them are only slightly different from the LO LFPT calculation, see the case with $\alpha=1/137$ in Fig. 3. Such an agreement confirms the equivalence of the tBLFQ approach and the LFPT calculation in the small coupling regime. Second, in the strong coupling regime, both the tBLFQ and the NLO LFPT calculations dramatically differ from the LO LFPT calculation, see the case with $\alpha=79/137$ in Fig. 3. Note that the period of the transition rate as a function of exposure time has changed due to higher order effects. The tBLFQ calculation should be regarded as a good approximation of the all order resummation results. Thus our comparison at strong coupling indicates that higher-order effects are significant for the interaction between the charged fermion and the electromagnetic field generated by an ultra-relativistic heavy ion. The plots in Fig. 3 also display the anticipated growth of higher-order effects with increasing atomic number of the nucleus.

IV Physical Observables

The electromagnetic field strength immediately after an ultra-relativistic heavy ion collision is proportional to the collision energy and reaches $m_{\pi}^{2}$ at RHIC and $10m_{\pi}^{2}$ at LHC. In addition, the field in the QGP medium could last up to a few fm/c. Such strong fields could lead to major modifications of physical observables, i.e.  the flow of the QGP, particle production, heavy quarkonium dissociation and so on, see Ref. Tuchin:2013ie for a review. In principle such modifications are within reach using the tBLFQ framework. However, as discussed in Sec. II.4, the fields generated by two colliding heavy ions moving in opposite directions have different dependence on light-front coordinates, thus it is a numerical challenge to study them simultaneously at this stage. In this paper, we study the real-time evolution of the momentum distribution of a fermion evolving inside the strong electromagnetic fields generated by one relativistic heavy ion as a demonstration of the tBLFQ formalism.

The momentum distributions of produced particles carry rich information about the collision process. In heavy ion collisions, the momentum distributions of various probes provide tomographic properties of the hot and dense medium created. Comparing to hadronic observables, electromagnetic probes, such as direct photons and dileptons, are valued for their greatly reduced final-state interactions. However, the strong magnetic field generated by the relativistic heavy ion Skokov:2009qp could significantly modify the configurations of charged particles, especially if they are produced immediately after the collision when the magnetic field is still strong Tuchin:2013ie . In this section we demonstrate real-time evolution of the momentum distribution for an electron evolving inside the strong electromagnetic field generated by an ultra-relativistic heavy nucleus using the tBLFQ formalism.

Consider the case in which a gold nucleus is moving along the positive $z$-axis with rapidity $y=5.3$. We perform our calculation starting at time $x^{\scriptscriptstyle+}=0$, at which the electron is in the BLFQ basis state $|\,\alpha_{0}\,\rangle$ with the following quantum number, $k=\frac{17}{2},n=0,m=0,\lambda=\frac{1}{2}$. As a result the initial longitudinal momentum of the electron is $p^{\scriptscriptstyle+}=17\pi/L$, and the initial transverse momentum distribution of the electron is,

[TABLE]

where $p^{\perp}\equiv|{\bm{p}}^{\perp}|$ and the transverse distribution is normalized as $\int\mathrm{d}p^{\perp}f_{0}(p^{\perp})=1$. The longitudinal length was chosen to be $L=10$ GeV*-1*, such that the longitudinal extent of the potential is limited. The HO parameter is set to be $b=1000m_{e}=0.511$ GeV in light of the fact that, in heavy ion collisions, the transverse momentum of the e.m. probes could reach a few GeV.

Using the tBLFQ formalism, we have access to the real-time evolution of the configuration of particles at the amplitude level. In this section, we illustrate this by showing the momentum probability distribution of the electron in both transverse and longitudinal directions as a function of time. We have considered the transition induced by both $\mathcal{W}_{1}$ and $\mathcal{W}_{2}$ terms of the Hamiltonian. We show snapshots of the transverse momentum probability distribution integrated over longitudinal momentum in Fig. 4 and the longitudinal momentum probability distribution in Fig. 6, simulated in tBLFQ and compared to the LO and NLO LFPT calculations using truncation parameters $N_{\text{max}}=32$ and $K=32$. We have checked that the momentum distribution is not sensitive to the truncation parameters $N_{\text{max}}$ and $K$. To be more specific, we observe that the calculation using $N_{\text{max}}=32$ and $K=32$ provides similar results comparing to calculation using $N_{\text{max}}=24$ and $K=24$ up to $x^{\scriptscriptstyle+}=50$ GeV*-1*. The numerical scheme adopted in the tBLFQ simulation is the MSD6 scheme, see Appendix B for details.

We show the transverse momentum distribution of the electron at $x^{\scriptscriptstyle+}=20$, $30$, $40$ and $50$ GeV*-1* from top to bottom in Fig. 4. The solid black curve is the initial transverse momentum distribution, which peaks around $0.7b$ with the peak value approximately equal to $1.7$. The tBLFQ simulation predicts that the transverse momentum distribution follows the LFPT predictions at first ($x^{\scriptscriptstyle+}\lesssim 20$ GeV*-1*). After being exposed to the intense field for a longer time, e.g., $x^{\scriptscriptstyle+}\gtrsim 30$ GeV*-1*, the peak value increases according to the tBLFQ simulation while both the LO and NLO LFPT predict that the value of the peak should decrease. Moreover, the tBLFQ simulation predicts that the position of the peak would be at a lower momentum comparing to the LFPT calculations, but with a smaller width. The different predictions by the tBLFQ simulation and the LFPT calculation could potentially be used as a quantitative observable for the higher order effects in this process. In Fig. 5, we show a snapshot of the transverse momentum distribution for $b<p_{\perp}<4b$ on a semi-log scale for the same process in Fig. 4 at $x^{\scriptscriptstyle+}=50$ GeV*-1*. It shows that the probability to find the electron with larger transverse momentum is significantly higher comparing to the initial distribution, if the electron has been exposed to the intense field for a sufficient amount of time. It is apparent because the electron has been excited to higher radial states from the initial $n=0$ state.

We show the longitudinal momentum distribution of the electron at $x^{\scriptscriptstyle+}=20$, $30$, $40$ and $50$ GeV*-1* from top to bottom in Fig. 6. Initially, the longitudinal momentum of the electron is $p^{\scriptscriptstyle+}=17\pi/L$. Since the strong e.m. field generated by the heavy ion has a non-trivial longitudinal momentum distribution, the evolving amplitude for the electron receives contributions showing that it is both accelerated as well as decelerated in the longitudinal direction. We discuss results of the tBLFQ simulation first. After being exposed to the strong field for $20$ GeV*-1*, the probability to find the electron in the initial longitudinal configuration has decreased to about $88\%$, and it is most likely to be found in states with longitudinal momentum adjacent to the initial momentum, e.g., the probabilities to find the electron with $p^{\scriptscriptstyle+}=19\pi/L$ and $p^{\scriptscriptstyle+}=15\pi/L$ are about $7\%$ and $1.5\%$, respectively. The transition rate to the $p^{\scriptscriptstyle+}=19\pi/L$ state is larger, because the kinetic energy difference between the $p^{\scriptscriptstyle+}=17\pi/L$ state and the $p^{\scriptscriptstyle+}=19\pi/L$ state is smaller comparing to the kinetic energy difference between the $p^{\scriptscriptstyle+}=17\pi/L$ state and the $p^{\scriptscriptstyle+}=15\pi/L$ state for the same transverse momentum distribution, owing to the fact that the kinetic energy of the electron is inversely proportional to the longitudinal momentum. The probability to find the electron in the initial longitudinal configuration continues decreasing with increasing exposure time in the intense field. At $x^{\scriptscriptstyle+}=50$ GeV*-1*, the probability to find the electron in the initial longitudinal configuration has decreased to about $75\%$. The probabilities to find the electron in states with $p^{\scriptscriptstyle+}=13\pi/L$ and $p^{\scriptscriptstyle+}=21\pi/L$ are about $1\%$ and $2\%$, respectively, which are not negligible. The probabilities to find the electron with $p^{\scriptscriptstyle+}=19\pi/L$ and $p^{\scriptscriptstyle+}=15\pi/L$ have increased to $14\%$ and $4\%$, respectively. The probabilities to find the electron in other longitudinal momentum states are also building up over time.

Comparing to the tBLFQ simulation, the LO LFPT calculation underestimates the depletion of the initial longitudinal momentum state. At $x^{\scriptscriptstyle+}=50$ GeV*-1*, the probabilities to find the electron in the initial longitudinal configuration, states with $p^{\scriptscriptstyle+}=19\pi/L$ and $p^{\scriptscriptstyle+}=15\pi/L$ are about $83\%$, $6\%$ and $6\%$, respectively. On the other hand, the NLO LFPT calculation predicts slightly higher transition probability to other longitudinal momentum states. At $x^{\scriptscriptstyle+}=50$ GeV*-1*, the probabilities to find the electron in the initial longitudinal state, and the states with $p^{\scriptscriptstyle+}=19\pi/L$ and $p^{\scriptscriptstyle+}=15\pi/L$ are about $71\%$, $16\%$ and $4\%$, respectively.

V Conclusions and Outlook

In this paper we use the recently formulated time-dependent basis light-front quantization formalism to study the influence of an electromagnetic field generated by an ultra-relativistic nucleus on a charged fermion field. We show that the information of the system is accessible at any intermediate time at the amplitude level using the tBLFQ approach. We show the BLFQ basis calculation is compatible with the momentum basis calculation with sufficiently large truncation parameters. Further, we calculate the transition of an electron influenced by the field generated by an ultra-relativistic nucleus and it shows agreement with light-front perturbation theory when the atomic number of the ultra-relativistic nucleus is small. We find that higher-order contributions are significant for nuclei with a large atomic number. We then demonstrate that the real-time evolution of the momentum distribution of an electron evolving inside the strong electromagnetic field can be calculated non-perturbatively using the tBLFQ approach.

Next, we plan to apply tBLFQ to QCD processes in high energy nuclear collisions, following the same procedures we presented in this investigation. For example, we can take the semi-analytic solution of the quasi-classical early time gluon field created in high energy nuclear collisions as background fields Fries:2006pv ; Chen:2015wia , and study the evolution of quarks and gluons in this field. Thus we could calculate high energy jet and heavy quark modification by the early time gluon field which could lead to improved understanding of jet quenching and heavy quark physics in heavy ion collisions Qin:2015srf ; Rapp:2008qc . We also plan to apply the tBLFQ formalism to electron-ion collisions Accardi:2012qut , to study the evolution of quarks and gluons in classical color fields. For example, we could study diffractive processes in electron-ion collisions using the dipole picture along with the classical description of a high energy nucleus from CGC, such that the dipole cross section and vector meson light-front wavefunction are obtained in a unified formalism Chen:2016dlk . The advantages of the tBLFQ framework are distinctive: it is a non-perturbative, first-principles numerical scheme; the calculation is at the amplitude level thereby incorporating quantum interference effects; and we can naturally extend our calculation to higher Fock-sectors as well as go beyond the Eikonal approximation.

Further improvement of tBLFQ relies on the development of BLFQ itself. For example, progress on implementing a sector-dependent renormalization scheme within the BLFQ framework Zhao:2014hpa will allow inclusion of higher Fock-sectors in our calculation; and what is more important, a proper renormalization scheme enables us to study various processes involving particle production and annihilation. For instance, we would be enabled to apply tBLFQ to the lepton pair and quarkonium production processes in ultra-peripheral heavy ion collision Baltz:2007kq ; Afanasiev:2009hy ; Abbas:2013oua . As computing technology advances, we envision that tBLFQ formalism will become a tool with increasing utility.

Acknowledgements.

We thank P. Maris for valuable remarks and constructive criticisms. We acknowledge A. Ilderton for fruitful discussions and M. Li for checking some equations in the paper. This work was supported by the Department of Energy under Grant Nos. DE-FG02-87ER40371 and DESC0008485 (SciDAC-3/NUCLEI). X. Zhao is supported by the new faculty startup funding by the Institute of Modern Physics, Chinese Academy of Sciences. We acknowledge computational resources provided by the National Energy Research Scientific Computing Center (NERSC), which is supported by the Office of Science of the U.S. Department of Energy under Contract No. DE-AC02-05CH11231.

Appendix A Conventions

The conventions we use in this paper are summarized in this section. Light-front coordinates are related to covariant coordinates $(x^{0},x^{1},x^{2},x^{3})$ as follows,

[TABLE]

with $x^{\scriptscriptstyle+}$ regarded as light-front time, and $x^{\scriptscriptstyle-}$ is the longitudinal coordinate. $x^{\perp}=(x^{1},x^{2})$ are the transverse coordinates. Non-vanishing elements of the metric tensor are $g^{+-}=g^{-+}=2$ and $g^{11}=g^{22}=-1$.

The basis states in the transverse direction are the eigenstates of the following two-dimensional harmonic oscillator (2D-HO) Hamiltonian

[TABLE]

in which $M$ and $\Omega$ are the mass and frequency of the oscillator. The characteristic scale of the 2D-HO depends only on a combination of these two parameters which we denote as $b=\sqrt{M\Omega}$. The eigenstates of Eq. (27) have two quantum numbers, the radial excitation $n$, and angular momentum projection $m$. The eigenenergy of a state with quantum number $n$ and $m$ is $E_{n,m}=(2n+|m|+1)\Omega$.

The basis wavefunctions in polar coordinates $(\rho,\phi)$, with $x^{1}=\rho\cos\phi$ and $x^{2}=\rho\sin\phi$, are,

[TABLE]

where the radial part $f_{nm}(\rho)$ is expressed by generalized Laguerre polynomials, $L^{|m|}_{n}(b^{2}\rho^{2})$ as,

[TABLE]

and angular part is

[TABLE]

A Fourier transform of the HO coordinate space wavefunctions immediately gives HO wavefunctions in momentum space,

[TABLE]

with

[TABLE]

and

[TABLE]

The coordinate and momentum space wavefunctions Eqs. (28) and (31) differ only in an overall coefficient if expressed as dimensionless parameter $b\rho$ and $p/b$.

The mode expansion of the fermion field operators in the BLFQ basis is,

[TABLE]

where ${\sf p\cdot x}=\frac{1}{2}p^{+}x^{\scriptscriptstyle-}-\bm{p}^{\perp}\cdot\bm{x}^{\perp}$ is the 3-product for the spatial components of $p^{\mu}$ and $x^{\mu}$. The creation operators $b^{\dagger}_{\bar{\alpha}}$ and $d^{\dagger}_{\bar{\alpha}}$ create electrons and positrons respectively with quantum numbers $\bar{\alpha}=\{k,n,m,\lambda\}$. They obey the anti-commutation relations

[TABLE]

We use the following (chiral) spinor representation, with helicity $\lambda=\pm 1/2$,

[TABLE]

Appendix B Multistep Differencing Scheme

Various schemes have been proposed for solving Eq. (4) numerically. For example, the Crank-Nicholson scheme (CN), which is unconditionally stable and accurate up to $(P^{\scriptscriptstyle-}\delta x^{\scriptscriptstyle+})^{2}$. However it is an implicit scheme which requires matrix inversion, a non desirable feature demanding tremendous computation efforts. There is also the Chebyshev scheme which approximates the exponential function by a Chebyshev polynomial expansion, which is stable and accurate however the intermediate wavefunctions are not available.

The multistep differencing scheme is an extension of the Euler scheme, which is stable and accurate while providing intermediate wavefunctions. The second order differencing scheme (MSD2) Askar:1978 relates the state at $x^{+}{+}\delta x^{+}$ to those at $x^{+}$ and $x^{+}{-}\delta x^{+}$ via

[TABLE]

It is conditionally stable if $|P^{\scriptscriptstyle-}_{\text{max}}|\delta x^{+}<1$, where $P^{\scriptscriptstyle-}_{\text{max}}$ is the largest (by magnitude) eigenvalue of $P^{\scriptscriptstyle-}$ when $P^{\scriptscriptstyle-}$ is time-independent Iitaka:1994 .

It has been shown that higher order multistep differencing scheme can provide much higher accuracy with some increase of computation efforts Iitaka:1994 . The fourth order scheme MSD4,

[TABLE]

which is accurate up to $(P^{\scriptscriptstyle-}\delta x^{+})^{4}$ and stable if $|P^{\scriptscriptstyle-}_{\text{max}}|\delta x^{+}<0.4$. The sixth order scheme MSD6,

[TABLE]

it is accurate up to $(P^{\scriptscriptstyle-}\delta x^{+})^{6}$ and stable when $|P^{\scriptscriptstyle-}_{\text{max}}|\delta x^{+}<0.1$.

The accuracy of MSD6 scheme for the calculation performed in Section III can be checked by comparing the evolution of eigenstates $|\,B\,\rangle$ of the Hamiltonian in Eq. (11), i.e., $P^{\scriptscriptstyle-}|\,B\,\rangle=P^{\scriptscriptstyle-}_{B}|\,B\,\rangle$, which is just a phase factor $\exp(iP^{\scriptscriptstyle-}_{B}\Delta x^{\scriptscriptstyle+})$. Note we use $|\,B\,\rangle$ to avoid confusion with eigenstate $|\,\beta\,\rangle$ of $P^{\scriptscriptstyle-}_{\text{QED}}$. The calculations in the BLFQ basis are tested to be accurate up to $4$ significant figures at $x^{\scriptscriptstyle+}=50$ GeV*-1* by successively halving the time increment.

Appendix C QED Hamiltonian in the BLFQ basis

The Hamiltonian relevant to the calculation we perform is summarized in Eq. (11). They are fermion kinetic energy $T_{f}$, vertex interaction between (anti-) fermion and background fields $\mathcal{W}_{1}$, and instantaneous-fermion interaction between (anti-)fermion and background fields $\mathcal{W}_{2}$. $T_{f}$ in BLFQ basis has been discussed in Zhao:2013cma . Here we outline how to express $\mathcal{W}_{1}$ and $\mathcal{W}_{2}$ algebraically in BLFQ basis following the convention in Appendix A.

The vertex interaction between fermion and background fields $\mathcal{A}^{\mu}$ as in Eqs. (14) is,

[TABLE]

where $\bar{\alpha}_{1},\bar{\alpha}_{2}$ are the quantum numbers associated with the field operators $\Psi$ and $\bar{\Psi}$ respectively. $k=(k^{\scriptscriptstyle+},\bm{k}^{\perp})$ is the momentum 3-vector of the background fields. We have in the transverse direction the 2D Dirac delta function and the Kronecker delta for the discretized longitudinal momentum. The instantaneous-fermion interaction between fermion and background fields $\mathcal{A}^{\mu}$ as in Eqs. (14) is,

[TABLE]

and we have used the HO wavefunction representation of the Dirac delta function in the transverse direction,

[TABLE]

so the integration with respect to $\bm{p}^{\perp}_{1},\bm{p}^{\perp}_{2},\bm{k}_{1}^{\perp},\bm{k}_{2}^{\perp}$ can be factorized. The delta function is exact only if $n_{i}$ and $m_{i}$ are summed over all possible values. However in LFPT, $\mathcal{W}_{2}$ has a singularity and it is canceled by the second order vertex interaction. Which implies we should regard $n_{i}$ and $m_{i}$ as quantum numbers of an intermediate fermion line and should be subjected to the same $N_{\text{max}}$ truncation constraint.

We list the spinor background potential contraction for different fermion helicities in Table 1. Note that the exponential phase factors in Eq. (14) are suppressed in the table.

Integration over the product of more than one highly oscillatory, 2D-HO wavefunctions, as in Eq. (40,41) can be simplified by applying the Talmi-Moshinsky transformation to the 2D-HO wavefunctions Davies . Eventually we are dealing with integration,

[TABLE]

which can be calculated as a finite-term summation using a series expansion of the Laguerre polynomials.