Lattice QED in an external magnetic field: Evidence for dynamical chiral symmetry breaking

We simulate QED in a strong constant homogeneous external magnetic field on a euclidean space-time lattice using the Rational Hybrid Monte Carlo method, developed for simulating lattice QCD. Our primary goal is to measure the chiral condensate in the limit when the input electron mass $m$ is zero. We observe a non-zero value, indicating that the external magnetic field catalyzes chiral symmetry breaking as predicted by approximate truncated Schwinger-Dyson methods. Such behaviour is associated with dominance by the lowest Landau level which causes the effective dimensional reduction from $3+1$~dimensions to $1+1$ dimensions for charged particles (electrons and positrons) where the attractive forces of QED can produce chiral symmetry breaking with a dynamical electron mass and associated chiral condensate. Since our lattice simulations use bare (lattice) parameters, while the Schwinger-Dyson analyses work with renormalized quantities, direct numerical comparison will require renormalization of our lattice results.


I. INTRODUCTION
Theoretical studies of free electrons in external electromagnetic fields were some of the earliest applications of relativistic quantum mechanics, [1][2][3], and exhibited some of the features of quantum electrodynamics, such as the Sauter-Schwinger effect [1,3] (the instability of the vacuum in a strong electric field, with the production of electron-positron pairs).[For a review of those field configurations where an exact closed form solution as simple integral or series is known, see, for example, the review article by Gerald Dunn.[4]] QED (quantum electrodynamics) in strong external electromagnetic fields is of relevance to high and medium energy physics, laser and accelerator physics, astrophysics and cosmology, and condensed matter physics.See, for example, V. Yakimenko et.al. [5], for a discussion with examples of the physics associated with QED in strong electromagnetic fields.
Current and future experiments at SLAC, LBNL and ELI, which collide the light from petawatt lasers with electron/positron beams or plasmas, produce environments with electromagnetic fields strong enough to produce quantum and QED effects [6][7][8].This requires electromagnetic fields which approach or exceed the critical values E cr = m 2 /e or B cr = m 2 /e where m and e are the electron mass and charge respectively.Another source of such strong fields is in the charged particle beams at future accelerators and their interactions, where fields could even be strong enough that the electron loop expansion (the last remnant of perturbation theory) breaks down [5].Certain compact astronomical objects, in particular those which are sources of X-and gamma-ray emissions and are identified as neutron stars are believed to have very strong surface magnetic fields (magnetars).See for example [9].Some have postulated the presence of very strong magnetic fields in the early universe as the source of the magnetic fields observed in the universe today.
We simulate lattice QED in external electromagnetic fields using the methods developed for simulating lattice QCD.Although some of the more interesting physics from an experimental point of view such as the Sauter-Schwinger effect -the production of electronpositron pairs from the unstable vacuum -are produced by strong external electric fields, this electric field makes the action complex, the imaginary part describing vacuum decay.
Hence standard lattice simulations cannot be applied.We therefore start with simulations of QED in external magnetic fields where the action is real and bounded below.We simulate using a non-compact gauge action and staggered fermions, using the RHMC algorithm to allow tuning to a single electron 'flavour'.For details see section 2.
We start with considering the case of 'free' electrons in a constant (in space and time) magnetic field with magnitude B, comparing our lattice results with the known exact solutions [3,4] to determine the range of B or more precisely eB (e is the electron charge) over which there is good agreement.
Next we perform simulations of full lattice QED at a near-physical electron charge α = e 2 /4π = 1/137, on a 36 4 lattice with safe (36 m >> 1, m << 1) electron masses m = 0.1, m = 0.2, comparing observables with those for free (α = 0) electrons in a magnetic field.We then perform simulations at B = 0 for a range of α on 36 4 lattices with m = 0.1, checking that the gauge action per site's α dependence is consistent with perturbation theory.
One of the theoretically most interesting predictions for QED in an external magnetic field is that the presence of this external field breaks chiral symmetry at m = 0.This manifests itself by giving a small dynamical mass to the electron which in turn gives rise to a non-zero chiral condensate.These predictions were obtained using truncated Schwinger-Dyson analyses, where the effects of the truncations required to obtain results, are difficult to estimate.The dynamical mass production was predicted in [10][11][12][13][14][15][16][17][18][19][20].Estimates of the chiral condensates are given in [21,22] For a good review with a more complete set of references see [23].For a more recent review, see, for example [24].It is therefore important to check these predictions using methods whose errors are easier to estimate and which allow systematic improvements.Lattice QED simulations are of this nature.However, the best estimates of Gusynin, Miransky and Shovkovy [12] predict a dynamical electron mass at our chosen eB = 0.4848... and a near-physical α = 1/137, m dyn ∼ 3 × 10 −35 , far below anything we could measure on the lattice.Therefore we choose a stronger electron charge α = 1/5, where the predicted m dyn ≈ 3 × 10 −4 .Here our simulations show evidence of chiral-symmetry breaking at a level 1 to 2 orders of magnitude greater than the 'best' Schwinger-Dyson results, however, our lattice QED results are for bare quantities, whereas those using Schwinger-Dyson analyses are for renormalized quantities.
Preliminary results were presented at Lattice 2021, Lattice 2022 and Lattice 2023 [25][26][27] Section 2 defines Lattice QED in an external magnetic field.Section 3 compares the effects of an external magnetic field on electrons on the lattice with those in the continuum.Section 4 presents simulations of Lattice QED at α = 1/137 in an external magnetic field.Section 5 describes simulations of Lattice QED at α = 1/5 in an external magnetic field.Section 6 presents discussions and conclusions.

II. LATTICE QED IN AN EXTERNAL MAGNETIC FIELD
In this section we describe the lattice transcription of QED in an external magnetic field used for our simulations.We use the non-compact gauge action where n is summed over the lattice sites and µ and ν run from 1 to 4 subject to the restriction.
The functional integral to calculate the expectation value for an observable O(A) where M = M † M and M is the staggered fermion action in the presence of the dynamic photon field A and external photon field A ext describing the magnetic [16] field B (or rather eB).M is defined by where the operator D µ is defined by and η µ are the staggered phases.Note that this treatment of the gauge-field-fermion interactions is compact and so has period 2π in the gauge fields.
We implement the RHMC simulation method of Clark and Kennedy [28], using a rational approximation to M −1/8 and rational approximations to M ±1/16 .To account for the range of normal modes of the non-compact gauge action, we vary the trajectory lengths τ over the range, of the periods of the modes of this gauge action.Here N µ is the length of the lattice in the µ direction which is chosen to be that of maximum extent [29].(Note that we typically only change the trajectory length after a trajectory has been accepted.) A ext are chosen in the symmetric gauge as while A ext,3 (n) = A ext,4 (n) = 0 [30].In practice we subtract the average values of A ext , µ from these definitions.This choice produces a magnetic field eB in the +z direction on every 1, 2 plaquette except that with i = N 1 , j = N 2 , which has the magnetic field eB(1 Because of the compact nature of the interaction, requiring eBN 1 N 2 = 2πn for some integer n = 0, 1, ......N 1 N 2 /2 makes the value of this plaquette indistinguishable from eB.Hence At the end of each accepted trajectory, we subtract the multiple of 2π/N µ from each A µ which reduces the magnitude of the lattice average of said A µ to lie in the range Since this is a gauge transformation it does not change any physics.In addition we transform to Landau gauge.Both these operations aim to prevent the gauge fields from becoming too large.
One of the observables we calculate is the electron contribution to the effective action per For ln(M) we use a (30,30) rational approximation to the logarithm.Here we use the Chebyshev method of Kelisky and Rivlin [31].While this has worse errors than a Remez approach, it preserves some of the properties of the logarithm itself, and is applicable on the whole complex plane cut along the negative real axis.(This would be important if we had an electric field, in which case det M would be complex).

III. 'FREE' ELECTRONS IN AN EXTERNAL MAGNETIC FIELD
We restrict ourselves to considering an external magnetic field of strength B which is constant in space and time and oriented in the +z direction.Classically 'free' electrons in such a magnetic field traverse helical orbits.The motion parallel to the magnetic field is free, while that orthogonal to the magnetic field is circular and hence bound.
Quantum mechanics restricts the motion perpendicular to the field to a discrete set of transverse energy levels known as the Landau levels [32].As B increases, the radii of these orbits decreases.The radius of the lowest Landau level is 1/ √ eB.This leads to an effective dimensional reduction from 3+1 to 1+1 dimensions for charged particles in a large magnetic field.The energy of the n-th Landau level is: where n = 0, 1, 2....The degeneracy of each n = 0 energy level is eB 2π A xy and that of each n > 0 energy level is eB π A xy where A xy is the area in the xy plane.We now calculate the chiral condensate ψψ on the lattice for masses m = 0.1 and m = 0.2 on a 36 4 lattice, values which we use for our initial QED simulations at α = 1/137, over a range of allowed eB values.This is to test the range of applicability of the lattice approach against the known continuum results.The lattice chiral condensate where M is defined in section 2. Because A ext is independent of z and t, we only need to calculate the trace over one xy plane.In fact, for m = 0.1 and m = 0.2, all terms in the trace are almost identical.We compare this with the known continuum result: where ψψ | eB=0 is taken from the lattice.This is necessary because ψψ | eB=0 depends on the UV regulator, which is different on the lattice than in the continuum.(Note that we have checked that we get the same result from equation 8 and from summing over the normal modes on the lattice.) Figure 1 compares the chiral condensates for free electrons in an external magnetic field B for m = 0.1 and for m = 0.2 as functions of eB on a 36 4 lattice, and compares them with the known continuum results.We conclude that the lattice results are in acceptable agreement with the continuum results for eB < ∼ 0.63 for each mass.We also calculate the fermion effective action/site on the lattice which we compare with the known continuum result: Again we replace the divergent part of this quantity in the continuum version with that of the lattice version to take into account the difference between the continuum and lattice regulators.Since the quadratically divergent part of L f , L f | eB=0 is a constant (independent of eB) and the logarithmically divergent term is proportional to (eB) 2 , while the leading contribution of the integral (the finite part) is of order (eB) 4 , the separation of the divergent and the finite parts is straight forward in principle.
Figure 2 shows the fermion contributions to the effective action on a 36 4 lattice as functions of eB, comparing the lattice and continuum results for m = 0.1 and m = 0.2.Again, we see good agreement over a range of eB values at least as large as for the chiral condensate.We simulate lattice QED on a 36 4 lattice using the approach presented in section 2 at lattice (bare) coupling α = 1/137 and lattice (bare) masses m = 0.1 and m = 0.2.At this α, the difference between bare and renormalized parameters is small (at most a few percent), and will therefore be ignored.Hence we consider these simulations to be performed at the We store a gauge configuration every 100 trajectories for future analyses.Note that, since QED probably does not have a UV completion, the action we choose should be considered to define an effective field theory.Its form is chosen to generate results consistent with QED perturbation theory with a lattice regulator.
In figure 3 the chiral condensates as functions of eB for α = 1/137 are compared to those for free electrons in a magnetic field (α = 0).In both cases the condensate for α = 1/137 lies above that for free electrons.This was to be expected, since the attractive force between electrons and positrons in QED is predicted to enhance the chiral condensate.
A truncated Schwinger-Dyson approach indicates that chiral symmetry breaking with a dynamical electron mass proportional to √ eB and a non-zero chiral condensate proportional to (eB) 3/2 survives in the m → 0 limit, however small α might be.However, at α = 1/137, as stated in the introduction, the dynamical electron contribution to the electron mass over the range of eBs accessible to these simulations is negligable so that these measurements could be checked by lattice perturbation theory.
Figure 4 compares the effective fermion actions/site for QED with α = 1/137, to its free-electron value as functions of eB for m = 0.1 and m = 0.2.The effective actions for QED are well above those without QED.We are interested in finding evidence of chiral symmetry breaking for QED in (strong) external magnetic fields, which survives in the m → 0 limit.Since Schwinger-Dyson analyses predict that this manifests itself as a dynamical electron mass, which at α = 1/137 and our chosen magnetic field eB = 0.4848... is m dyn ≈ 3 × 10 −35 , in the m = 0 limit, far below anything we could measure on the lattice, we simulate at a far stronger bare coupling, α = 1/5.For a renormalized α = 1/5, the predicted m dyn ≈ 3 × 10 −4 [12] equation ( 51).
This leads to ψψ ≈ 1.2 × 10 −4 , [22] equation (B4), which should be measurable.Of course, since our α = 1/5 is the bare(lattice) α, there is no guarantee that it will give a measurable result, we only know this a postiori.
First, we need to determine that α = 1/5 is still in the range of α values where perturbation theory holds.We do this by performing simulations over a range of α values from α = 0 to α = 1/5.These simulations were performed at eB = 0 and m = 0.1 Our chosen observable for performing this test is the gauge Lagrangian L γ , the gauge action/lattice site.
The fact that L γ is well approximated by a low order polynomial in α with coefficients of decreasing magnitude over the range 0 ≤ α ≤ 0.2 is consistent with the claim that α = 0.2 lies in the perturbative domain for eB = 0. We have also performed simulations at α = 1, m=0.1 for which L γ (1) = 1.49047(3), that does not appear to be consistent with perturbation theory, although one would need more simulations in the range 0.2 ≤ α ≤ 1.0 to check this.
Having chosen α = 1/5 = 0.2 for our lattice QED simulations at strong coupling, we first perform simulations at eB = 0. Since we do not expect any chiral symmetry breaking breaking in the m → 0 limit, the chiral condensate at non-zero mass should be proportional to mΛ 2 probably with logarithmic corrections, where Λ is the momentum cutoff.On the lattice, since we use units where the lattice spacing is 1, Λ = ±π in each of the 4 directions x, y, z, t.This means that the trace of the propagator, which defines the chiral condensate, is expected to be dominated by large momenta and hence short distances, and should not be sensitive to the size of the lattice.Therefore we should be able to calculate the condensate on relatively small lattices, even in the limit as m → 0. We simulate over a range of electron masses 0.001 ≤ m ≤ 0.2 on a 36 4 lattice, noting that for the lowest mass m = 0.001, mN µ = 0.036 << 1, far outside the range of 'safe' values for which would require mN µ >> 1. ments of chiral condensates on the lattice for free fermions in an external magnetic field show appreciable departures from known continuum values.Here we need to make measurements at m values small enough that our 36 4 lattice is too small to yield infinite lattice values for the chiral condensate.However, the assumption that for large eB only the lowest Landau levels (LLL) make significant contributions to physics means that as long as the lattice projection in the xy plane is considerably larger than that of the LLL, whose radii are ≈ 1/ √ eB, then the chiral condensates should not depend on N x = N y , no matter how small m becomes.For this reason we fix N x = N y = 36 or N x = N y = 18 for our simulations.
On the other hand, we expect that, at small mass, the chiral condensates will depend on the lattice extents in the z and t directions.
We perform simulations at a selection of m values in the range 0.001 ≤ m ≤ 0.2, starting with a 36 4 lattice at each mass.At the lower masses we then increase the lattice size to 36 2 × 64 2 .Because at m = 0.025, the increase in the chiral condensate in going from a 36 4 lattice to a 36 2 × 64 2 lattice is very small, we conclude that a 36 4 lattice would have been adequate, and that it is unnecessary to increase the lattice sizes for any m > 0.025.For m = 0.0125 there is a small but significant increase in the chiral condensate in going from a 36 4 lattice to a 36 2 × 64 2 lattice, from which we conclude that a 36 4 lattice is too small, but a 36 2 × 64 2 lattice is probably adequate.For m = 0.005 we conclude that a 36 2 × 96 2 lattice is adequate while for m = 0.001 we needed the 18 2 × 128 2 lattice.
Figure 7 shows the chiral condensates as functions of mass at α = 1/5 for eB = 2π × 100/36 2 = 2π × 25/18 2 = 0.4848....This strongly suggests that the condensates approach a finite, non-zero limit as m → 0. From the obvious curvature of this graph (through the uppermost point at each m), a straight line through the points for the 2 lightest will pass through m = 0 at a point marking an upper bound to the m = 0 condensate.This value is 3.977... × 10 −3 .To estimate a lower bound, we need to fit a smooth curve to the condensates for the smallest masses for which we have simulations.This requires choosing a functional form which displays increasing curvature as m → 0, and fitting it to the condensates for the lowest masses.The first form we choose is: Of course, a direct comparison is not possible because the parameters in the lattice are bare (lattice) parameters, while those in the Schwinger-Dyson estimates are renormalized quantities.We discuss this further in the Discussion and Conclusions section.

VI. DISCUSSION AND CONCLUSIONS
We simulate lattice QED in a strong external magnetic field using the RHMC method developed for simulating lattice QCD.Much of our effort is aimed at seeking evidence for chiral symmetry breaking in the limit that the lattice (bare) mass approaches zero, catalyzed by the external magnetic field, as predicted by less reliable truncated Schwinger-Dyson analyses.Such chiral symmetry breaking is facilitated in part by the electrons and positrons preferentially occupying the lowest Landau level (LLL) with its radius 1/ √ eB due to the external magnetic field causing an effective dimensional reduction from 3 + 1 to 1 + 1 for charged particles.The attractive force due to QED then produces chiral symmetry breaking.
The broken chiral symmetry manifests itself by producing a dynamical mass for the electron and a non-zero chiral condensate.
Note that the parameters and fields in the action are bare lattice quantities.This compares with the renormalized condensate of ≈ 1.2 × 10 −4 from the best Schwinger-Dyson estimates.
Since the lowest non-trivial contributions to the renormalization condensates evaluated at p 2 = m 2 dyn are typically O(α log(Λ 2 /m 2 dyn )) ≈ 1-3 where Λ ≈ π is the momentum cutoff on the lattice, renormalization cannot be ignored.This contrasts with our simulations with α = 1/137, m = 0.1, 0.2, where the renormalization constants are evaluated at p 2 = m 2 , so the lowest non-trivial order corrections are typically O(α log(Λ 2 /m 2 ) ≈ 0.05, 0.04) and can be ignored if one can tolerate errors of a few percent (in actual fact for small external magnetic fields the errors are somewhat smaller than this).There are two ways that one might check the consistency between the Schwinger-Dyson approach and the Lattice QED simulations.The first is to follow the methods developed for renormalizing lattice QCD.The dimensional reduction from 3 + 1 to 1 + 1 dimensions caused by the external magnetic field adds extra complications.The second approach is to repeat the Schwinger-Dyson analysis on the lattice action based on methods developed for lattice perturbation theory, in terms of bare parameters.Of course, it would be best to use both methods.
Experience with lattice QCD indicates that we should expect large taste breaking effects from the use of staggered fermions (and, in particular, rooted staggered fermions).See for example the review article [33], and its guide to the literature.In lattice QCD these effects are, to lowest order, O(a 2 Λ 2 QCD ) where a is the lattice spacing.For lattice QED in an external magnetic field these effects should be, to lowest order O(a 2 eB) or, since we have chosen a = 1, O(eB).Hence we should be able to reduce taste breaking by reducing the external magnetic field.We are therefore repeating our simulations with a smaller magnetic field.Another way one can see that this should decrease errors is to note that with a smaller magnetic field, the lowest Landau level covers more lattice sites so that discretization errors should be reduced.With the smaller magnetic field, we should be able to test that the chiral condensate scales as (eB) 3/2 as predicted.
It would be interesting to apply lattice simulations to QED with N f = 2 flavours in an external magnetic field, where such chiral symmetry breaking would imply spontaneous dynamical breaking of flavour symmetry with associated Goldstone bosons.Note that these Goldstone bosons are uncharged and hence reside in 3 + 1 dimensions and therefore spontaneous breaking of continuous chiral symmetry is allowed despite the electrons being restricted to 1 + 1 dimensions.
One of the first calculations planned for stored configurations is to measure the effect that QED in an external magnetic field has on the coulomb field of a point charge in said magnetic field.It has been predicted that the coulomb field will be partially screened and distorted by the presence of the external magnetic field [34][35][36][37].This effect can and will be measured on our stored configurations at α = 1/137 where renormalization can be safely ignored, with expected errors of at most a few percent.Such an analysis will be performed using Wilson loops, possibly with smearing.
We plan to calculate electron propagators on our stored configurations with α = 1/5, eB = 0.4848... on lattice sizes from 36 2 × 64 2 to 18 2 × 128 to extract the electron mass extrapolated to m = 0.When configurations from our current and ongoing simulations at eB = 0.1163... become available we will then be able to check that m dyn ∝ √ eB, and extract the electron wave-function renormalization constant.We will also attempt to calculate the photon propagator, although this will be considerably more difficult, since its falloff with separation of its endpoints will be statistics limited.
The one parameter left to calculate will be the renormalized coupling constant which, unless this can be calculated perturbatively, might well require additional simulations.
Of at least as much interest is the behaviour of QED in an external electric field, which exhibits the Sauter-Schwinger effect.As pointed out earlier, this has a complex action and traditional lattice simulations which rely on importance sampling cannot be used.A good starting point is to simulate lattice QED with both external electric and magnetic fields where these fields are obtained by a boost from a pure external magnetic field, since the physics should be the same as that for this magnetic field alone.One method we might try is a Complex Langevin (CLE) simulation.This has some chance of succeeding for QED in external electric fields or electric and magnetic fields, whereas it failed for QCD in a quarknumber chemical potential, because for pure non-compact U(1) lattice gauge theory, which is a free field theory and hence a collection of harmonic oscillators, the real trajectory is an attractive fixed-point of the CLE, while for pure compact SU(3) lattice gauge theory it is a repulsive fixed point.We will first need to check if the real trajectory remains an attractive fixed point for the U(1) lattice gauge theory CLE simulation when the electron fields are included making it lattice QED, but without the external electric fields.If so, we plan to try CLE when the external electric or electric and magnetic fields are added.

FIG. 1 :
FIG. 1: Electron chiral condensates ψψ as functions of eB, comparing the continuum and lattice results for a) m = 0.1 and b) m = 0.2.

FIG. 2 :
FIG. 2: -Electron effective actions/site −L f as functions of eB, comparing the continuum and lattice results for a) m = 0.1 and b) m = 0.2.
physical electron charge.For m = 0.1 we simulate over a range of allowed eB values in the interval 0 ≤ eB ≤ 2π × 160/36 2 = 0.7757..., while for m = 0.2 we use a selection of allowed eB values in the range 0 ≤ eB ≤ 2π × 200 = 0.9696....At each eB we run for 12500 trajectories of lengths randomly chosen from the periods of the modes of the gauge action.

FIG. 5 :
FIG. 5: Photon Lagrangian density as a function of α for eB = 0.The curve is the fit

Figure 6 FIG. 6 :
Figure6shows the chiral condensate ψψ as a function of m on a 364 lattice for α = 1/5.This graph also shows the value of the chiral condensate at the lowest mass, m = 0.001, from a simulation on a larger (48 4 ) lattice, which shows no sign of any appreciable dependence on lattice size.In fact, ψψ (m = 0.001) = 4.3259(7) × 10 −4 on a 364 lattice compared with ψψ (m = 0.001) = 4.3294(7) × 10 −4 on a 48 4 lattice.A simple linear fit to the lowest 2 masses (m = 0.001,m = 0.005), which because of the curvature of this graph should yield an upper bound to the value at m = 0, predicts ψψ m=0 = 5(1) × 10 −7 , a value which curvature can easily lower to zero.Hence the chiral condensate at m = 0 for α = 1/5 and eB = 0 is consistent with zero.To search for evidence of chiral symmetry breaking at m = 0 catalyzed by a strong magnetic field we perform simulations at α = 1/5 and eB = 2π × 100/36 2 = 2π × 25/18 2 = 0.4848..., relatively large, while being significantly below eB = 0.63 above which measure-