Scaling dimensions in QED$_3$ from the $\epsilon$-expansion

We study the fixed point that controls the IR dynamics of QED in $d = 4 - 2\epsilon$. We derive the scaling dimensions of four-fermion and bilinear operators beyond leading order in $\epsilon$-expansion. For the four-fermion operators, this requires the computation of a two-loop mixing that was not known before. We then extrapolate these scaling dimensions to $d = 3$ to estimate their value at the IR fixed point of QED$_3$ as function of the number of fermions $N_f$. The next-to-leading order result for the four-fermion operators corrects significantly the leading one. Our best estimate at this order indicates that they do not cross marginality for any value of $N_f$, which would imply that they cannot trigger a departure from the conformal phase. For the scaling dimensions of bilinear operators, we observe better convergence as we increase the order. In particular, $\epsilon$-expansion provides a convincing estimate for the dimension of the flavor-singlet scalar in the full range of $N_f$.


Introduction
antum Electrodynamics (QED) in 3d is an asymptotically-free gauge theory, which becomes strongly interacting in the IR. When the U (1) gauge eld is coupled to an even number, 2N f , of complex two-component fermions, and the Chern-Simons level is zero, the theory is parity invariant and has an SU (2N f ) × U (1) global symmetry. For large N f the theory ows in the IR to an interacting conformal eld theory (CFT) that enjoys the same parity and global symmetry.
e CFT observables are then amenable to perturbation theory in 1/N f ; this has been done for scaling dimensions [1][2][3][4][5][6][7][8][9][10][11], two-point functions of conserved currents [12][13][14], and the free energy [15]. e IR xed point is expected to persist beyond this large-N f regime, but not much is known about it. Ref. [16] employed the conformal bootstrap approach to derive bounds on the scaling dimensions of some monopole operators. Another method to study the small-N f CFT is the -expansion, which exploits the existence of a xed point of Wilson-Fisher-type [17] in QED continued to d = 4 − 2 dimensions. When 1 we can access observables via a perturbative expansion in and subsequently a empt an extrapolation to = 1 2 . e -expansion of QED was employed to estimate some scaling dimensions [18,19], the free energy F [20], and the coe cients C T and C J [14]. In particular, ref. [18] considered operators made out of gauge-invariant products of either four or two fermion elds.
Four-fermion operators are interesting because of the dynamical role they can play in the transition from the conformal to a symmetry-breaking phase, which is conjectured to exist if N f is smaller than a certain critical number N c f [21][22][23][24]. In fact, the operators with the lowest UV dimension that are singlet under the symmetries of the theory are four-fermion operators.
If for small N f they are dangerously irrelevant, i.e., their anomalous dimension is large enough for them to ow to relevant operators in the IR, they may trigger the aforementioned transition [7,25,26]. 1 e one-loop result of ref. [18] led to the estimate N c f ≤ 2. Bilinear operators, i.e., operators with two fermion elds, are interesting because they are presumably among the operators with lowest dimension. For instance, when continued to d = 3, the two-form operators Ψγ [µ γ ν] Ψ become the additional conserved currents of the SU (2N f ) symmetry, of which only a SU (N f ) subgroup is visible in d = 4 − 2 . is leads to the conjecture that their scaling dimension should approach the value ∆ = 2 as → 1 2 , which was tested at the one-loop level in ref. [18].
In order to assess the reliability of the -expansion in QED, and improve the estimates from the one-loop extrapolations, it is desirable to extend the calculation of these anomalous dimensions beyond leading order in . is is the purpose of the present paper. Let us describe the computations we perform and the signi cance of the results.
We rst consider four-fermion operators. In the UV theory in d = 4 − 2 , there are two such operators that upon continuation to d = 3 match with the singlets of the SU (2N f ) symmetry. We compute their anomalous dimension matrix (ADM) at two-loop level by renormalizing oshell, amputated Green's functions of elementary elds with a single operator insertion. As we discuss in detail in a companion paper [30], knowing this two-by-two ADM is not su cient to obtain the O( 2 ) scaling dimensions at the IR xed point. We also need to take into account the full one-loop mixing with a family of in nitely many operators that have the same dimension in the free theory. ese operators are of the form where n is an odd integer, and Γ n µ 1 ...µn ≡ γ [µ 1 . . . γ µn] is an antisymmetrized product of gamma matrices. All the operators in this family except for the rst two, i.e., n = 1, 3, vanish for the integer values d = 4 and d = 3, but are non-trivial for intermediate values 3 < d < 4. For this reason they are called evanescent operators. Taking properly into account the contribution of the evanescent operators, via the approach described in ref. [30], we obtain the next-toleading order (NLO) scaling dimension of the rst two operators. We then extrapolate to = 1 2 using a Padé approximant, leading to the result presented in subsection 5.2 and summarized in gure 2. e deviation from the leading order (LO) scaling dimension is considerable for small N f , indicating that at this order we cannot yet obtain a precise estimate for this observable of the three-dimensional CFT. Taking, however, the NLO result at face value, we would conclude that the four-fermion operators are never dangerously irrelevant. is resonates with recent results that suggest that QED 3 is conformal in the IR for any value of N f . Namely, refs. [31][32][33] argued based on 3d bosonization dualities [34][35][36][37], that for N f = 1 the SU (2) × U (1) symmetry is in fact enhanced to O(4) (this is related to the self-duality present in this theory [38]). Also, a recent la ice study [39] found no evidence for a symmetry-breaking condensate (for previous la ice studies see refs. [40][41][42]). We then consider the bilinear "tensor-current" operators of the form ΨΓ n µ 1 ...µn Ψ , (1.2) for n = 0, 1, 2, 3. We obtain their IR scaling dimension up to O( 3 ) using the three-loop computations from ref. [43]. Having these higher-order results, we are in the position to employ di erent Padé approximants to estimate errors and test the convergence as we increase the order. As mentioned above, in the limit d → 3 the operators with n = 1, 2 approach conserved currents of the SU (2N f ) symmetry. Indeed, we show in subsection 5.3 (see gure 4) that the extrapolated scaling dimension of the two-form operators approaches the value ∆ = 2 as we increase the order. As d → 3, the operators with n = 0, 3 approach scalar bilinears, which are either in the adjoint representation of SU (2N f ) or are singlets. For the singlet scalar, which is continued by a bilinear with n = 3, the results of various extrapolations we perform are all close to each other (see gure 5), indicating that -expansion provides a good estimate for this scaling dimension in the full range of N f . For the adjoint scalar, di erent components are continued by operators with either n = 0 or n = 3, giving two independent extrapolations at each order in . As expected, we nd that the two independent extrapolations approach each other as we increase the order (see gure 5). e rest of the paper is organized as follows: in section 2 we set up our notation and describe the xed point of QED in d = 4 − 2 ; in section 3 we present the computation of the two-loop ADM of the four-fermion operators, and then the result for their scaling dimension at the IR xed point in d = 4 − 2 ; in section 4 we present the same result for the bilinear operators; in section 5 we extrapolate the scaling dimensions to d = 3, and plot the resulting dimensions as a function of N f for the various operators we consider; nally in section 6 we present our conclusions and discuss possible future directions. In the appendices we collect additional material and some useful intermediate results.

QED in d = 4 − 2
We consider QED with N f Dirac fermions Ψ a , a = 1, . . . , N f , of charge 1. e Lagrangian is with the covariant derivative de ned as Summation over repeated avor indices is implicit. We work in the R ξ -gauge, de ned by adding the gauge-xing term We collect the Feynman rules in appendix A. e algebra of the gamma matrices is {γ µ , γ ν } = 2η µν , with η µν η νρ = δ ρ µ and δ µ µ = d. We will employ some useful results on d-dimensional Cli ord algebras from ref. [44]. We normalize the traces by Tr[1] = 4, for any d. For d = 3, Ψ a decomposes as giving 2N f complex two-component 3d fermions ψ i , i = 1, . . . , 2N f , all with charge 1. Correspondingly, the gamma matrices decompose as µ } µ=1,2,3 are two-by-two 3d gamma matrices. In d = 4, the global symmetry preserved by the gauge-coupling is In d = 4 − 2 , evanescent operators violate the conservation of the nonsinglet axial currents [45], so only the diagonal subgroup SU (N f ) is preserved. In d = 3, this symmetry enhances to SU (2N f ) × U (1).
We de ne α ≡ e 2 16π 2 and denote bare quantities with a subscript "0". e renormalized coupling is given by where the renormalization constant Z α (α, ) absorbs the poles at = 0, and µ is the renormalization scale. e beta function reads In Minimal Subtraction (MS), β depends only on α and not on . e MS QED β function is known up to four-loop order for generic N f [46,47] Using eqs. (2.7) and (2.9) we nd that in d = 4 − 2 the theory has a xed point at with ζ(n) the Riemann zeta function. Our convention for renormalizing elds is By the Ward Identity, Z A = Z −1 α . For our computations we need the eld-renormalization of the fermion up to two-loop order. In MS and generic R ξ -gauge it reads (2.12)

Operator mixing
To compute the anomalous dimension of local operators O i , we add these operators to the Lagrangian and compute their renormalized couplings C i at linear level in the bare ones (2.14) Z j i are the mixing renormalization constants from which we obtain the ADM Like β, γ does not depend on in the MS scheme. We introduce the following notation for the coe cients of the expansion in α and e most direct way to compute the mixing Z j i is to renormalize amputated one-particleirreducible Green's functions with zero-momentum operator insertions and elementary elds as external legs. Alternatively, one can renormalize the two-point functions of the composite operators. e former method has two main advantages. e rst is that to extract n-loop poles only n-loop diagrams need to be computed. e second is that we can insert the operators with zero-momentum.
is makes higher-loop computations more tractable. e disadvantage is that o -shell Green's functions with elementary elds as external legs are not gauge-invariant, so some results in the intermediate steps of the calculation are ξ-dependent, which is why we need to include the ξ-dependent wave-function renormalization of external fermions. In addition, operators that vanish under the equations of motion (EOM) enter the renormalization of such o -shell Green's functions. We refer to the la er as EOM-vanishing operators.
In the next section, we consider composite operators given by scalar quadrilinear and bilinear operators in the fermion elds. We rst present the computation of the two-loop anomalous dimension of the four-fermion operators and use it to obtain the O( 2 ) IR scaling dimension at the xed point. Next, we employ the already existing results of the three-loop anomalous dimension of bilinear operators [43] to obtain their IR dimension to O( 3 ).

Four-fermion operators in d = 4 − 2
In this section, we present the computation of the ADM of the four-fermion operators at the two-loop level. Here and in the following Γ n µ 1 ...µn ≡ γ [µ 1 . . . γ µn] , with the square brackets denoting antisymmetrization, which includes the conventional normalization factor 1 n! . In d = 4, the operators in eq. (3.1) are the only two operators with scaling dimension 6 at the free xed point that are singlets under the global symmetry SU (N f ) L × SU (N f ) R . We focus on these avor-singlet operators because, as explained in the introduction, we are interested in understanding whether or not they are relevant at the IR xed point. e calculation of the ADM for avor-nonsinglet operators is actually simpler because it involves a subset of the diagrams. We report the result for some nonsinglet operators in appendix C.
In d = 4 − 2 , insertions of Q 1 and Q 3 in loop diagrams generate additional structures that are linearly independent to the Feynman rules of Q 1 and Q 3 . To renormalize the divergences proportional to such structures, we need to enlarge the operator basis. It is most convenient to de ne the complete basis by adding operators that vanish for → 0, and hence are called evanescent operators, as opposed to Q 1 and Q 3 that we refer to as physical operators. ere is an in nite set of such evanescent operators. One choice of basis for them is with n an odd integer ≥ 5. e terms proportional to the arbitrary constants a n and b n are of the form × a physical operator; they parametrize di erent possible choices for the basis of evanescent operators.
For the computation of the ADM we adopt the subtraction scheme introduced in refs. [48,49]. Since this is the most commonly used scheme for applications in avor physics, we refer to it as the avor scheme. We label indices of the ADM using odd integers n ≤ 1, so that n = 1, 3 correspond to the physical operators, eq. (3.1), and n ≥ 5 to the evanescent operators, eq. (3.2). e ADM up to two-loop order is 2 where  Notice that the invariant (Q 1 , Q 3 ) block of γ (2,0) depends on the coe cients a 5 and b 5 , which parametrize our choice of basis. is dependence can be understood as a sign of schemedependence [51]. Clearly, this implies that the scaling dimensions at O( 2 ) are not simply obtained from the eigenvalues of this invariant block, as also its eigenvalues depend on a 5 and b 5 . e additional contribution that cancels this basis-dependence originates from the O( ) term γ (1,−1) in the one-loop ADM. Such O( ) terms are indeed induced in every scheme that contains nite renormalizations, such as the avor scheme. For a thorough discussion of the scheme/basis-dependence and its cancellation we refer to ref. [30].
ere are a few non-trivial ways of partially testing the correctness of the two-loop results: i) We performed all computations in general R ξ gauge. is allowed us to explicitly check that the mixing of gauge-invariant operators indeed does not depend on ξ.
ii) All the two-loop counterterms are local, i.e., the local counterterms from one-loop diagrams subtract all terms proportional to 1 log µ in two-loop diagrams.
iii) e 1 2 poles of the two-loop mixing constants satisfy the relation where β (1,0) is the one-loop coe cient of the beta-function. is is equivalent to the -independence of the anomalous dimension [52].
In the next two subsections, we discuss the renormalization of the one-and two-loop Green's functions from which we extract the relevant entries of the mixing matrix Z -and ultimately the ADM entries in eqs. (3.4), (3.5), and (3.6)-and some technical aspects of the two-loop computation. A reader more interested in the results for the scaling dimensions may proceed directly to section 3.4.

Operator basis
As argued in section 2.1, in general we need to consider also EOM-vanishing operators when renormalizing o -shell Green's functions. Moreover, in our computation we adopt an IR regulator that breaks gauge-invariance, so we also need to take into account some gauge-variant operators. Below we list all operators that, together with (Q 1 , Q 3 ) and {E n } n≥5 , enter the renormalization of the two-loop Green's functions we consider:

EOM-vanishing operators
ere is a single EOM-vanishing operator, N 1 , that a ects the ADM at the one-loop level and another one, N 2 , that a ects it at the two-loop level. ey read Additionally, there are EOM-vanishing operators that are only necessary to close the basis of independent Lorentz structures for certain Green's functions. For completeness, we list them here Here / D ≡ γ µ D µ and the arrow indicates on which eld the derivative is acting.
Gauge-variant operators Renormalization constants subtract UV poles of Green's functions. It is thus essential to ensure that no IR poles are mistakenly included in the renormalization constants. In practice, this means that an energy scale must be present in dimensionally regularized integrals. Otherwise, UV and IR contributions cancel each other and the result of the loop integral is zero in dimensional regularization [45].
One possibility to introduce a scale is to keep the external momentum in the loop integral. However, i) such loop integrals are more involved than integrals obtained by expanding in Table 1: A summary of the Green's functions we consider. e loop order (L-loop) refers to the α L contribution to the corresponding Green's function (second column). e third column contains the mixing renormalization constants that the given Green's function depends on. e last column contains the ones we extract in each case. powers of external momenta over loop momenta, and ii) keeping external momenta does not necessarily cure all the IR divergences, e.g., diagrams with gluonic snails in non-abelian gauge theories. Another possibility for QED would be to introduce a mass for the Dirac fermions. e drawback in this case is that we would have to consider many more EOM-vanishing operators.
Instead, we apply the method of "Infrared Rearrangement" [53,54]. is method consists in rewriting the massless propagators as a sum of a term with a reduced degree of divergence and a term depending on an arti cial mass, m IRA . Section 3.3 contains some more details about the method. e caveat is that the method violates gauge invariance in intermediate steps of the computation. All breaking of gauge invariance is proportional to m 2 IRA and explicitly cancels in physical quantities. However, to restore gauge-invariance, also gauge-variant operators proportional to m 2 IRA need to be consistently included in the computation. Fortunately, due to the factor of m 2 IRA , at each dimension there are only a few of them. At the dimension-four level, there is a single operator generated, i.e., the photon-mass operator: At the dimension-six level, there are more operators, but only one, P, enters our ADM computation because Q 1 and Q 3 mix into it at one-loop. It reads

Renormalizing Green's functions
In this subsection, we highlight the relevant aspects in the computation of the renormalization constants Z j i , from which we extracted the ADM presented above, via the renormalization of amputated one-particle irreducible Green's functions.
For each Green's function we need to specify the operator we insert and the elementary elds on the external legs. In our case, the external legs are either four elementary fermions, or two fermions and a photon, or two photons. At the tree-level, a Wick contraction with the elementary elds de nes a vertex structure for each operator. We denote the ΨΨΨΨ structures with S, the ΨΨA µ ones withS, and the A µ A ν one withŜ. An additional subscript indicates the operator associated to a given structure. e representation in terms of Feynman diagrams is We collect all structures that enter the computation in appendix A.
In what follows, we refer to as a sum over a speci c subset of Feynman diagrams: i) All these diagrams have a single insertion of the operator O. ii) ey are dressed with interactions such that they contribute at O(α L ). In particular, we include all counterterm diagrams proportional to eld and charge renormalization constants, but we do not include diagrams that contain mixing constants. We keep those separate to demonstrate how we extract them. iii) e subscript S indicates that out of this sum of diagrams we only take the part proportional to the structure S. In short, the notation of eq. (3.14) denotes the L-loop insertion of O projected on S, including contributions from eld and charge renormalization constants.
As an illustration of the notation we show in gure 1 a small subset of the Feynman diagrams for the non-trivial case of N 1 (2) S , withS any of the structures in eq. (A.10). Notice that since N 1 is a linear combination of terms with di erent elds, see eq. (3.8), its eld and charge renormalizations depend on the part we insert, namely Next we derive the conditions on the Green's functions that determine the mixing constants. For transparency we frame the constant(s) that we extract from a given condition. In table 1 we summarize which Green's functions we consider, on which mixing renormalization constants they depend, and which one we extract in each case. For brevity we use the following shorthand notation: We collect the results for the renormalization constants in appendix B. A µ A ν at one-loop At one-loop there is no insertion of any four-fermion operator that contributes to the Green's function with only two external photons. us Contrarily, one-loop insertions of four-fermion operator contribute to the ΨΨA µ Green's function. By expanding the diagram in the basis ofS structures, we determine the mixing into operators with a tree-level projection onto ΨΨA µ , namely N 1 and P. For the physical operators the conditions are In the rst line we use that Z Notice that in this case the mixing constants subtract nite terms, as required by the avor scheme we adopt.

ΨΨΨΨ at one-loop
Next, we compute the one-loop insertions in the ΨΨΨΨ Green's function. Firstly, we insert physical operators, i.e., Q, with the only non-vanishing N 1 QN 1 , which we have previously determined via the ΨΨA µ Green's function. Next, we insert evanescent operators. Again, the only di erence here is that their mixing constants into physical operators subtract nite pieces is completes the computation of all one-loop constants required to determine the mixing of physical operators at the two-loop level. Next, we renormalize the same Green's functions at the two-loop level.
At the two-loop order Q 1 and Q 3 insertions do contribute to the A µ A ν Green's function. ey can thus mix into the operator N 2 . Even though N 2 itself does not have a tree-level projection on physical operators, we need this mixing to extract the two-loop mixing of Q 1 and Q 3 into N 1 in the next step. e projection onto theŜ structure results in the condition

(3.23)
ΨΨA µ at two-loop Next we renormalize the ΨΨA µ Green's function at the two-loop level. We only need the two-loop mixing of physical operators into N 1 , because only N 1 has a tree-level projection onto Q 1 . To unambiguously determine the projection on the structureS N 1 , we have to x a basis of linear independent structures, which correspond to linearly independent operators. At this loop order, we nd that apart from N 1 we also need to include the operators N 3 and N 4 to project all generated structures. is projection is the only point in which these operators enter our computation. e niteness of the two-loop ΨΨA µ Green's function determines the two-loop mixing of physical operators into N 1 via 3 ΨΨΨΨ at two-loop Finally, we have collected all results necessary to renormalize the two-loop ΨΨΨΨ Green's function. e renormalization conditions for the mixing in the physical sector read We see here explicitly that, because N 1 has a tree-level projection onto Q 1 , we need Z

Evaluation of Feynman diagrams
Already at the two-loop level the number of Feynman diagrams entering the Green's functions is quite large. e present computation is thus performed in an automated setup. Firstly, the program QGRAF [55] generates all diagrams creating a symbolic output for each diagram. is output is converted to the algebraic structure of a loop diagram and subsequently computed using self-wri en routines in FORM [56]. e methods for the computation and extraction of the UV poles of two-loop diagrams are not novel and also widely used throughout the literature. Here, we shall only sketch the steps and mention parts speci c to our computation.
One major simpli cation of the computation comes from the fact that we can always expand the integrand in powers of external momenta over loop-momenta and drop terms beyond the order we are interested in. For instance, for the ΨΨΨΨ Green's function all external momenta can be directly set to zero, while for the ΨΨA µ one we need to keep the external momenta up to second order to obtain the mixing into N 1 (seeS N 1 in eq. (A.10)).
A er the expansion, all propagators are massless so the resulting loop-integrals vanish in dimensional regularization. To regularize the IR poles and perform the expansion in external momenta we implement the "Infrared Rearrangement" (IRA) procedure introduced in refs. [53,3 Note that N1 (1) , as N1 has two Feynman rules. 54]. In IRA, an -in our case massless-propagator is replaced using the identity where p is the loop momentum, q is a linear combination of external momenta, and m IRA is an arti cial, unphysical mass. We see that the rst term in the decomposition contains the scale m IRA and carries no dependence on external momenta in its denominator. In the second term, the original propagator reappears, but thanks to the additional factor the overall degree of divergence of the diagram is reduced by one. When we apply the decomposition multiple times, we obtain a sum of terms with only loop-momenta and m IRA in the denominators plus a term proportional to 1 (p+q) 2 . is last term, however, can be made to have an arbitrary small degree of divergence. erefore, in a given diagram we can always perform the decomposition as many times as necessary until terms proportional to 1 (p+q) 2 are nite and can thus be dropped if we are interested in UV poles.
When applying IRA on photon propagators, the resulting coe cients of the poles are not gauge-invariant, because we drop the nite terms in the expansion of propagators. is is why some gauge-variant operators/counterterms enter in intermediate stages of the computation, for instance the operator P. Such operators are always proportional to m 2 IRA and so only a small number of them enters at each dimension. For more details on the prescription we refer to the original work [54].
e IRA procedure results in integrals with denominators that i) are independent from external momenta, and ii) contain the arti cial mass m IRA . We can always reduce these integrals to scalar "vacuum" diagrams by contracting them with metric tensors and solving the resulting system of linear equations, e.g., see ref. [54]. is tensor reduction reduces all integrals to oneand two-loop scalar integrals of the form with the integers n 1 , n 2 , n 3 ≥ 1, and m 1 = 0. e one-loop integral can be directly evaluated, whereas all two-loop integrals can be reduced to a few master integrals using the recursion relation in ref. [57]. In fact, in our case m 1 = m 2 = m IRA and the use of recursion relations is not required.
In the evaluation of the Feynman diagrams, we use the Cli ord algebra in d dimensions for i) the evaluation of traces with gamma matrices when the diagram in question has closed fermion loops, and ii) the reduction of the Dirac structures to the operator structures S orS listed in appendix A.

Anomalous dimensions at the fixed point
where not required otherwise . (3.31) Note that the physical-physical block is not invariant at order 2 , because there are non-zero entries (γ * ) n1 and (γ * ) n3 for all n ≥ 5.
We are interested in nding the rst two eigenvalues of γ * up to order 2 . ey determine the scaling dimensions of the corresponding eigenoperators at the IR xed point. We denote these scaling dimensions by with i = 1, 2 and ∆ UV ( ) = 6−4 . To compute the rst two eigenvalues we have truncated the problem to include a large but nite number of evanescent operators. Taking a su ciently large truncation, the scheme/basis-dependence of the approximated result can be made negligible at the level of precision we are interested in (for details see ref. [30]). In table 2, we list the values of (∆ 1 ) i and (∆ 2 ) i for N f = 1, . . . , 10 a er we included enough evanescent operators such that the three signi cant digits listed remain unchanged. e table is the main result of this section. In section 5, we will use these results as a starting point to extrapolate the scaling dimensions to d = 3.

Bilinear operators in d = − 2
In this section we consider operators that are bilinear in the fermionic elds. e most generic bilinear operators without derivatives are  Table 2: e values of the one-loop (∆ 1 ) i and the two-loop (∆ 2 ) i coe cients de ned in eq. In d = 4−2 , the conservation of the nonsinglet axial currents is violated by evanescent operators [45], and thus only the diagonal SU (N f ) is a symmetry. On the other hand, the CFT in d = 3 is expected to enjoy the full SU (N f ) L × SU (N f ) R symmetry, which is actually enhanced to SU (2N f ) × U (1). erefore, in continuing the operators of eq. (4.1) to d = 3, we nd that the ones with γ 5 are in the same multiplets of the avor symmetry as those without. So even though their scaling dimensions can di er as a function of , the enhanced symmetry entails that they should agree when = 1 2 . Since the operators with γ 5 do not provide new information about the 3d CFT, and the 't Hoo -Veltman prescription makes computations technically more involved, we restrict our discussion here to operators without γ 5 . As a future direction, it would be interesting to test this prediction of the enhanced symmetry by comparing the scaling dimensions of operators with γ 5 a er extrapolating to d = 3 at su ciently high order. We also restrict the discussion to operators with n ≤ 3, because the others are evanescent in d = 3.
e anomalous dimension of bilinear operators without γ 5 has been computed for a generic gauge group at three-loop accuracy in ref. [43]. For our U (1) gauge theory we substitute C A = 0 and C F = T F = 1. Moreover, there is a di erence in the normalization convention for the anomalous dimension, so that γ here = 2γ there . Under SU (N f ) each operator decomposes into a singlet and an adjoint component, respectively. A priori, the two components can have di erent anomalous dimensions. e di erence between the singlet and the adjoint originates from diagrams in which the operator is inserted in a closed fermion loop. When the operator has an even number of gamma matrices, the closed loop gives a trace with an odd total number of gamma matrices, which vanishes. So for even n there is no di erence between the singlet and the adjoint, i.e., they have the same anomalous dimension.
Below we collect the results for n ≤ 3.  Two-form:

Scalar:
Three-form: In d = 4 these three-form operators are Hodge-dual to axial currents. Actually, the fact that they do not get an anomalous dimension at one-loop, as seen from the equations above, is related to this. However, Hodge-duality cannot be de ned in d = 4 − 2 and the anomalous dimensions start to di er from those of the axial current at the two-loop level.
is exhausts the list of bilinears without γ 5 that ow to physical operators as d → 3. In section 5.3 we discuss which operators of the CFT in d = 3 are continued by the operators above, and extrapolate the above results to obtain estimates for their scaling dimensions.

Padé approximants
A computation of a certain order in provides an approximation to the observable, e.g. the scaling dimension ∆, in terms of a polynomial  Padé (1,1) 6.86 6.52 6.35 6.25 6.19 6.15 6.12 6.10 6.08 6.07 Taking → 1 2 in this polynomial gives the " xed order" d = 3 prediction of the -expansion. Typically, the xed-order results show poor convergence as the order is increased. A standard resummation technique adopted for these kind of extrapolations is to replace the polynomial with a Padé approximant. e Padé approximant of order (k,l) is de ned as e coe cients c i and d i are determined by matching the expansion of eq. (5.2) with eq. (5.1). k + l must equal the order at which we are computing. Another condition comes from the fact that we are interested in the result for → 1 2 . In order for the -expansion to smoothly interpolate from = 0 to = 1 2 , an employable Padé approximant should not have poles for ∈ [0, 1 2 ] for the values of N f that we consider. In what follows, we show the predictions from a Padé approximation only if it does not contain any pole on the positive axis of for any value of N f = 1, . . . , 10.

Four-fermion operators as d → 3
In d = 3, the two four-fermion operators in the UV can be rewri en as where i = 1, . . . , 2N f . In this rewriting we see explicitly that these operators are singlets of SU (2N f ). We now evaluate the scaling dimensions (∆) 1 and (∆) 2 of the two corresponding IR eigenoperators, at NLO. For the NLO prediction we employ the Padé approximation of order (1,1). We list the values of the LO and NLO Padé (1,1) predictions for the values of N f = 1, . . . , 10 in table 3.
We visualize the results in gure 2. e dashed lines are the result of the one-loopexpansion computation. Indeed, as discussed in ref. [18], the one-loop approximation predicts that the lowest eigenvalue becomes relevant for N f < 3. e two-loop computation presented here changes this prediction. e two solid lines represent the NLO Padé (1,1) approximation to the two scaling dimensions. We observe that for no value of N f does the lowest eigenvalue reach marginality. We also see that the corrections to the LO result are signi cant, especially for small N f , i.e., N f = 1, 2. is means that for such small values of N f , NLO accuracy is not su cient to obtain a precise estimate for this scaling dimension. Nevertheless, at face value, the result of the two-loop -expansion suggests that QED 3 is conformal in the IR for any value of N f . Next, we comment on the relation of our result to the 1/N f -expansion in d = 3. At large µ ψ i , is set to zero by the EOM of the gauge eld, hence the operator Q 1 is an EOM-vanishing operator. However, besides Q 3 , there still is another avor-singlet scalar operator of dimension 4 for N f = ∞, namely F 2 µν . Q 3 and F 2 µν mix at order 1/N f [11]. Looking at the -expansion result in gure 2 we see that indeed only the lowest eigenvalue (∆) 1 (black lines) approaches 4 for large N f . e other scaling dimension (red lines) approaches 6 as N f → ∞, implying that the two eigenoperators cannot mix at large N f . is is consistent precisely because there is only one non-trivial singlet four-fermion operator at large N f . Its mixing with F 2 µν cannot be captured within the -expansion, because the UV dimension of F 2 µν di ers from that of a four-fermion operator in d = 4 − 2 . We can, however, test whether for any value of ∈ [0, 1 2 ] the lowest eigenvalue (∆) 1 , which starts o larger at = 0, crosses the dimension of F 2 µν . Such a level-crossing would require to revisit the extrapolation to = 1 2 and possibly a ect the estimate. e scaling dimension of F 2 µν in -expansion is with α * given in eq. (2.10) up to O( 4 ). At three-and four-loop order the only Padé approximation without poles in the positive real axis of is the order (2,1) and (2,2), respectively. In gure 3 we plot (∆) 1,2 and ∆(F 2 ) as a function of d for the representative cases of N f = 1, 2, fermion operators (black and red lines) and F 2 µν (blue lines) as a function of the dimension d, i.e., for ∈ [0, 1 2 ]. e le , center, and right panel show the result for the representative cases of N f = 1, 2, and 10, respectively. We observe that the N 3 LO Padé (2,2) prediction of ∆(F 2 ) never crosses the NLO Padé (1,1) prediction of (∆) 1 in the extrapolation region. and 10. We observe that the only case in which (∆) 1 crosses ∆(F 2 ) before d = 3 is when N f = 1 and when we employ N 2 LO Padé (2,1) to predict ∆(F 2 ). e N 3 LO Padé (2,2) prediction for N f = 1 does not cross (∆) 1 and the same holds for larger values of N f . erefore, at least at this order, F 2 µν should not play a signi cant role in obtaining the four-fermion scaling dimension.

Bilinears as d → 3
Next we consider bilinear operators in d = 3. In the UV, restricting to the ones without derivatives, the possibilities are Scalar: e subscript refers to the representation of SU (2N f ). e singlet is parity-odd. We can combine parity with an element of the Cartan of SU (2N f ), in such a way that one component of the adjoint scalar is parity-even. Since parity squares to the identity, this Cartan element can only have +1 and −1 along the diagonal, which up to permutations we can take to be the rst N f , and the second N f diagonal entries, respectively. With this choice, the parity-even bilinear is N f a=1 (ψ a ψ a − ψ a+N f ψ a+N f ). is is the candidate to be the "chiral condensate" in QED 3 [22]. Vector: e singlet is the current of the gauged U (1). When the interaction is turned on, it recombines with the eld strength and does not ow to any primary operator of the IR CFT.
e adjoint is the current that generates the SU (2N f ) global symmetry. erefore, we expect it to remain conserved along the RG and ow to a conserved current of dimension ∆ = 2 in the IR.
We now identify which d = 4 − 2 bilinears from section 4 approach the d = 3 bilinears above. Substituting the decomposition of eqs. (2.4) and (2.5), and also using 3d Hodge duality, we nd that In gure 4 we plot the extrapolations for the scaling dimension of the conserved avornonsinglet current B (1) adj as a function of N f . We observe that both N 2 LO Padé approximants are closer to 2 than the LO and NLO ones, and they remain close to 2 even for small values of N f . We consider this to be a successful test of the -expansion, which supports its viability as a tool to study QED 3 .     sing we nd good convergence behaviour between the NLO Padé (1,1) and the two N 2 LO Padé approximations. erefore, for this observable we are able to provide a rather convincing estimate. We do stress, however, that the comparison of the various approximations does not provide rigorous error estimates, since the error due to the extrapolation is not under control. For B (0) adj we have two di erent operators that provide a continuation to d = 4 − 2 . It is encouraging that as the order increases, the two resulting estimates approach each other. Even so, we nd that for small N f the N 2 LO Padé approximations are spread, so the -expansion at this order does not provide a de nite prediction. As N f increases the situation improves, namely all NLO and N 2 LO approximations begin to converge.
In table 4 we list the numerical values for the various estimates of the bilinear scaling dimensions for N f = 1, . . . , 10.
Next, we compare to the large-N f predictions for the scaling dimensions of the bilinears. e Padé approximants used to estimate the dimensions of B and compare the coe cient c (k,l) with its exact value obtained from the large-N f expansion,   is suggests that the extrapolation of the three-form may provide a be er estimate for the scaling dimension of the adjoint scalar at this order.

Conclusions and future directions
We employed the -expansion to compute scaling dimensions of four-fermion and bilinear operators at the IR xed point of QED in d = 4 − 2 . We estimated the corresponding value for the physically interesting case of d = 3. e results seem to con rm the expectations from the enhancement of the global symmetry as d → 3 (see gures 4 and 5). erefore, going beyond the leading order gave us more con dence that the continuation is sensible. At the same time, it appears that -with the exception of the scalar-singlet bilinear-to obtain precise estimates for the scaling dimensions for small values of N f requires even higher-order computations and perhaps more sophisticated resummation techniques (see for instance chapter 16 of ref. [60] and references therein). e computation of such higher orders in via the standard techniques used in the present work would require hard Feynman-diagram calculations.
On a di erent note, ref. [84] recently argued that QCD 3 with massless quarks undergoes a transition from a conformal IR phase, which exists for su ciently large number of avors, to a symmetry-breaking phase when N f ≤ N c f . is is analogous to the long-standing conjecture for QED 3 , and so four-fermion operators may play the same role.
erefore, at least for the case of zero Chern-Simons level, -expansion can be employed in a similar manner to estimate N c f . A LO estimate appeared in ref. [85]. In light of our results for QED 3 , it would be worth studying how this estimate is modi ed at NLO.
Acknowledgements: we thank Joachim Brod, Martin Gorbahn, John Gracey, Igor Klebanov, Zohar Komargodski, and David Stone for their interest and the many helpful discussions. We are also indebted to the Weizmann Institute of Science, in which this research began. Research at Perimeter Institute is supported by the Government of Canada through Industry Canada and by the Province of Ontario through the Ministry of Research & Innovation.

A Feynman rules
From the QED Lagrangian in R ξ -gauge, we obtain the Feynman rules ere is one additional counterterm coupling that we need to specify. It is a relic of the procedure with which we regulate IR divergences (see section 3.3), which essentially breaks gauge invariance. For this reason to consistently renormalize Green's function we need to include a counterterm analogous to a mass for the photon, i.e., Only the one-loop value of δm 2 IRA enters our computations. It reads To nd the EOM-vanishing operators at the non-renormalizable level we apply the EOM of the fermion and photon. ey read For brevity we use the shorthand notation γ µ D µ ≡ / D and use an arrow to indicate the direction in which the derivative in / D acts, i.e. / D and / D ≡ / D. We consider the Lagrangian with additional couplings proportional to the operators introduced in section 3.1 To compute the Green's function we need the Feynman rules of the operators we insert, as well as all the structures that we need to project the amplitude. For instance, to renormalize the Green's function of ΨΨA µ with one-loop insertions of Q 1 we need not only the Feynman rule of Q 1 , but also the ΨΨA µ structure of all operators that Q 1 generates at one-loop.

B Renormalization constants
In this appendix we list the mixing-renormalization constants of four-fermion operators. First we list the constants we need to compute the ADM of avor-singlet four-fermion operators, which we discussed in the main text, and subsequently the constants entering the computation of the ADM of avor-nonsinglet four-fermion operators, which we discuss in appendix C.

B.1 Flavor-singlet four-fermion operators
and for the evanescent operators they are with n an odd integer ≥ 5. To compute these constants for generic n we used Cli ord-algebra identities from ref. [44]. As explained in section 3.2, in the computation of the mixing at two-loop level more operators enter. e only one-loop mixings entering the computation, apart from those above, is the mixing of the physical four-fermion operators into the EOM-vanishing operator N 2 , and the gauge-variant operator P. e former vanish, i.e., with Q = Q 1 , Q 3 . We do not list the corresponding constants for the evanescent operators because they do not enter the two-loop computation of the mixing of physical operators.
In table 1 we summarised on which renormalization constants the Green's functions we computed depend on. We see that to determine the two-loop mixing of the four-fermion operators we rst need to determine the two-loop mixing of the physical operators into the two EOM-vanishing operators N 1 and N 2 . e corresponding constants read Z (2,2) Finally, the two-loop mixing constants of the two physical operators read

B.2 Flavor-nonsinglet four-fermion operators
e renormalization of the Green's functions with insertions of avor-nonsinglet four-fermion operators is analogous to the one with avor-singlets but less involved. eir avor-o -diagonal structure forbids them to receive contributions from any EOM-vanishing or gauge-variant operator at two-loop order. erefore, in this case we only need the mixing constants within the physical and evanescent sectors.
As in the avor-singlet case, the one-loop mixing is directly related to the one-loop anomalous dimensions of eqs. (C.4) and (C.5) via with O, O any physical or evanescent avor-nonsinglet four-fermion operator; the one-loop anomalous dimensions above are given in appendix C. Finally, the two-loop mixing constants  Table 5: ree signi cant digits of the one-loop, (∆ 1 ) i , and the two-loop, (∆ 2 ) i , contributions to the scaling dimension of the avor-nonsinglet four-fermion operators for various cases of N f . To obtain the two-loop (∆ 2 ) i values we implemented the algorithm to include the e ect of evanescent operators [30].

C Flavor-nonsinglet four-fermion operators
In the main part of this work we investigated bilinear and avor-singlet four-fermion operators. ere exist also four-fermion operators that are not singlets under avor. e ones we consider in this appendix are spanned by the basis E n = T ac bd (Ψ a Γ n µ 1 ...µn Ψ b )(Ψ c Γ n µ 1 ...µn Ψ d ) + a n Q 1 + b n Q 3 , (C. 3) with T ac db = T ca bd and T ac ad = T ab bd = 0. e computation of their ADM at one-and two-loop order entails only a subset of the Feynman diagrams needed for avor-singlet case and is actually less involved as discussed in appendix B. In this appendix we present their ADM and their scaling dimensions at the IR xed point in d = 4 − 2 , and use this to estimate the corresponding d = 3 observables.
(C.6) e part of the one-loop result that does not depend on a n and b n was rst computed in ref. [48].  Table 6: LO and either NLO Padé (1,1) or xed-order NLO predictions for the scaling dimension of the two avor-nonsinglet four-fermion operators at d = 3 for various values of N f . Only three signi cant digits are being displayed.