Two-loop evolution equations for flavor-singlet light-ray operators

TL;DR

This paper derives two-loop evolution equations for flavor-singlet light-ray operators in QCD, confirming previous results and enhancing understanding of operator mixing and conformal symmetry in quantum field theory.

Contribution

It introduces a novel method leveraging conformal invariance in non-integer dimensions to compute two-loop evolution kernels for flavor-singlet operators.

Findings

Calculated two-loop evolution kernels for flavor-singlet operators.

Confirmed previous evolution equations for generalized parton distributions.

Demonstrated the utility of conformal symmetry in higher-order QCD calculations.

Abstract

QCD in non-integer $d=4-2\epsilon$ space-time dimensions enjoys conformal invariance at the special fine-tuned value of the coupling. Counterterms for composite operators in minimal subtraction schemes do not depend on $\epsilon$ by construction, and therefore the renormalization group equations for composite operators in physical (integer) dimensions inherit conformal symmetry. This observation can be used to restore the complete evolution kernels that take into account mixing with the operators containing total derivatives from their eigenvalues (anomalous dimensions). Using this approach we calculate the two-loop (NLO) evolution kernels for the leading twist flavor-singlet operators in the position space (light-ray operator) representation. As the main result of phenomenological relevance, in this way we are able to confirm the evolution equations of flavor-singlet generalized hadron parton distributions derived earlier by Belitsky and M\"uller using a different approach.

Equations

$${\mathcal{O}}_g(z_1,z_2)$$

$${\mathcal{O}}_q(z_1,z_2)$$

$$[z_1,z_2]=P\exp \left\{ig{z}_{12}{\int }_0^{1}du{A}_{+}(z_{21}^u)\right\},$$

$$z_{12}=z_1-z_2,z_{21}^u=\bar{u}{z}_2+u{z}_1,\bar{u}=1-u.$$

$${\mathcal{O}}_{\chi }(z_1,z_2)=\sum _{n=1}^{\infty }\sum _{k=n}^{\infty }{\Psi }_{nk}^{(\chi )}(z_1,z_2){\mathcal{O}}_{\chi ,nk}(0),$$

$${\mathcal{O}}_{q,nk}={\partial }_{+}^k\bar{q}{C}_n^{(\frac{3}{2})}\left(\frac{{\overleftarrow{D}}_{+}-{\overrightarrow{D}}_{+}}{{\overleftarrow{D}}_{+}+{\overrightarrow{D}}_{+}}\right)q,$$

$$[{\mathcal{O}}_{\chi }(z_1,z_2)]=\sum _{n,k}{\Psi }_{nk}^{(\chi )}(z_1,z_2)[{\mathcal{O}}_{\chi ,nk}(0)].$$

$$\displaystyle\left(\Big{(}M\partial_{M}+\beta(a)\partial_{a}\Big{)}\delta^{\chi\chi^{\prime}}+\gamma^{\chi\chi^{\prime}}_{nn^{\prime}}(a)\right)[\mathcal{O}_{\chi^{\prime},n^{\prime}k}]=0\,,$$

$$\displaystyle a=\frac{\alpha_{s}}{4\pi}\,,\qquad\beta(a)=M\frac{da}{dM}=-2a\big{(}\epsilon-\gamma_{g}\big{)}=-2a\big{(}\epsilon+a\beta_{0}+a^{2}\beta_{1}+\ldots\big{)}$$

$$\displaystyle\beta_{0}=\frac{11}{3}C_{A}-\frac{2}{3}n_{f}\,,\qquad\beta_{1}=\frac{2}{3}\big{[}17C_{A}^{2}-5n_{f}C_{A}-3n_{f}C_{F}\big{]}\,.$$

$$\displaystyle\Big{(}M\partial_{M}+\beta(a)\partial_{a}\Big{)}[\mathcal{O}_{\chi}](z_{1},z_{2})$$

$${\mathbb{H}}_{\chi {\chi }^{\prime }}(a)=a{\mathbb{H}}_{\chi {\chi }^{\prime }}^{(1)}+a^2{\mathbb{H}}_{\chi {\chi }^{\prime }}^{(2)}+\dots .$$

$${\mathbb{H}}^{(1)}\left(\begin{array}{c}{\mathcal{O}}_q(z_1,z_2)\\ {\mathcal{O}}_g(z_1,z_2)\end{array}\right)$$

$$-{\int }_{\alpha \beta }\left(\begin{array}{cc}4C_F& 4n_f{z}_{12}(\bar{\alpha }\bar{\beta }+3\alpha \beta )\\ 8C_F{z}_{12}^{-1}& 16C_A(\bar{\alpha }\bar{\beta }+2\alpha \beta )\end{array}\right)\left(\begin{array}{c}{\mathcal{O}}_q(z_{12}^{\alpha },z_{21}^{\beta })\\ {\mathcal{O}}_g(z_{12}^{\alpha },z_{21}^{\beta })\end{array}\right)$$

$$+\left(\begin{array}{cc}2C_F& 0\\ -4C_F{z}_{12}^{-1}& 12C_A-2{\beta }_0\end{array}\right)\left(\begin{array}{c}{\mathcal{O}}_q(z_1,z_2)\\ {\mathcal{O}}_g(z_1,z_2)\end{array}\right).$$

$${\int }_{\alpha }\equiv {\int }_0^{1}d\alpha ,$$

$$S_{\alpha }^{(0)}{\mathbb{H}}^{(1)}={\mathbb{H}}^{(1)}{S}_{\alpha }^{(0)},$$

$$S_{\chi ,-}^{(0)}$$

$$S_{\chi ,0}^{(0)}$$

$$S_{\chi ,+}^{(0)}$$

$$j_q=1,j_g=3/2.$$

$$[{\mathcal{O}}_{\chi ,nk}(0)]={\hat{\mathcal{O}}}_{\chi ,nk}+{\mathcal{B}}_{\chi ,nk}+{\mathcal{E}}_{\chi ,nk},$$

$${\hat{\mathcal{O}}}_{\chi }(z_1,z_2)=\sum _{nk}{\Psi }_{nk}^{(\chi )}(z_1,z_2){\hat{\mathcal{O}}}_{\chi ,nk}(0).$$

$${\delta }_{\alpha }{\hat{\mathcal{O}}}_{\chi }(z_1,z_2)\equiv \sum _{nk}{\Psi }_{nk}^{(\chi )}(z_1,z_2){\delta }_{\alpha }{\hat{\mathcal{O}}}_{\chi ,nk}(0),$$

$${\delta }_{\alpha }{\hat{\mathcal{O}}}_{\chi ,nk}(0)=\sum _{{\chi }^{\prime },n}(c_{\alpha }{)}_{n{n}^{\prime }}^{\chi {\chi }^{\prime }}{\hat{\mathcal{O}}}_{{\chi }^{\prime },n^{\prime },k+\alpha }(0)+\dots ,$$

$${\delta }_{\alpha }{\hat{\mathcal{O}}}_{\chi }(z_1,z_2)=\sum _{{\chi }^{\prime }}{S}_{\alpha ,\chi {\chi }^{\prime }}{\hat{\mathcal{O}}}_{{\chi }^{\prime }}(z_1,z_2).$$

$$S_{-}$$

$$S_0$$

$$S_{+}$$

$$\Delta (a_{\ast })=a_{\ast }{\Delta }^{(1)}+\mathcal{O}(a_{\ast }^2)$$

$${\Delta }^{(1)}\left(\begin{array}{c}{\mathcal{O}}_q(z_1,z_2)\\ {\mathcal{O}}_g(z_1,z_2)\end{array}\right)$$

$$[\mathbb{H}(a_{\ast }),S_{+}(a_{\ast })]=0,$$

$$[S_{+}^{(0)},{\mathbb{H}}^{(1)}]$$

$$[S_{+}^{(0)},{\mathbb{H}}^{(2)}]$$

Full text

1 Introduction

The exciting possibility to study the three-dimensional proton structure (“tomographic imaging”) at the next generation lepton-hadron and lepton-nucleus colliders, most notably the EIC Accardi:2012qut , poses formidable theoretical challenges. To meet these challenges precise predictions in quantum chromodynamics (QCD) are necessary and require pushing the corresponding computational tools to the highest possible accuracy. In practice, however, there exists still a considerable gap between precision and theoretical rigor applied in data analyses at the energy frontier, e.g., by the LHC community searching for New Physics beyond the Standard Model, and the QCD studies at medium energy machines, where the main goal is an improved understanding of the strong interactions.

This concerns, in particular, differences in the accuracy when considering the scale dependence of parton distributions. Whereas the next-to-next-to-leading order (NNLO), i.e., three-loop, analysis of parton distributions and fragmentation functions is becoming standard Accardi:2016ndt , the analysis of deeply-virtual Compton scattering and, e.g., hard exclusive vector meson electro-production is still often based on the leading-order (LO) evolution, despite the complete next-to-leading order (NLO) evolution kernels being available for a long time Belitsky:1999hf . The necessity to close this gap is becoming gradually accepted so that the task to verify the results of Ref. Belitsky:1999hf by an independent calculation and to develop the techniques to push such calculations to NNLO accuracy are high on the agenda.

On the technical level, the challenge is that considering off-forward matrix elements one has to take into account mixing with the operators containing total derivatives. The direct calculation of the relevant Feynman diagrams is difficult as they involve two different momenta, so that the method of choice has been to use the constraints from conformal symmetry that allow one to “save” one loop in comparison to the direct computation. This approach was pioneered by Dieter Müller Mueller:1991gd ; Mueller:1993hg ; Mueller:1997ak and further refined in Refs. Belitsky:1998gc ; Belitsky:1998vj ; Belitsky:1999hf . The existing results for the two-loop (NLO) evolution kernels Belitsky:1999hf have been obtained using this technique.

In Ref. Braun:2013tva we proposed a somewhat different implementation of the same idea, based on the exact conformal symmetry of QCD in $d=4-2\epsilon$ dimensions at the critical point. Another difference to Müller’s approach is our use of the position space (light-ray operator) representation that allows one to switch rather easily to momentum fraction space, avoiding the generally nontrivial problem of the restoration of the evolution kernels from the mixing matrices of local operators. Evolution kernels in position space are also interesting on their own in connection with lattice QCD calculations of Euclidean “observables” that can be factorized in terms of parton distributions, see e.g., Refs. Braun:2007wv ; Ji:2013dva ; Ma:2017pxb . Such calculations currently attract a lot of attention.

This modified approach was tested in  Braun:2013tva on several examples to two- and three-loop accuracy for scalar theories. In  Braun:2014vba ; Braun:2016qlg ; Braun:2017cih we used it to calculate the NLO (two-loop) Braun:2014vba and NNLO (three-loop) Braun:2016qlg ; Braun:2017cih evolution kernels for the leading-twist quark-antiquark flavor-nonsinglet operators. In this work we address the flavor-singlet sector. This case is technically more complicated and also new questions arise concerning the role of gauge-noninvariant operators in conformal Ward identities. These issues will be discussed in what follows. Our main result is the derivation of the complete set of NLO evolution kernels in the position space representation for the leading twist gluon and (C-even) quark-antiquark operators. Expanding these kernels at small field separations we reproduce the results for the mixing matrices for flavor-singlet local operators given in Ref. Belitsky:1998gc . We also demonstrate how the evolution kernels in momentum fraction space can be obtained from the position space expressions by simple integration. We have verified numerically that the resulting expressions agree with the kernels that are implemented in the NLO evolution code for generalized parton distributions (GPD) by Freund and McDermott Freund:2001rk ; Freund:2001hd which is based on the analytic expressions derived in Ref. Belitsky:1999hf .

2 Flavor-singlet light-ray operators

The general formalism of evolution equations for light-ray operators goes back to Ref. Balitsky:1987bk and is described in detail in  Braun:2014vba ; Braun:2016qlg ; Braun:2017cih so that we will be brief in what follows. In this work we study the evolution equations for the flavor-singlet $C$-even twist-two light-ray operators (LROs) which are defined as 111All equations are written assuming Euclidean space.

[TABLE]

Here $n^{\mu}$ is an auxiliary light-like vector, $n^{2}=0$, the “plus” projection stands for $F_{+\mu}=n^{\nu}F_{\nu\mu}$, and $z_{1},z_{2}$ are real numbers. We will often omit the factor $n^{\mu}$ in the argument specifying the field position and use a short-hand notation $q(z)\equiv q(zn)$, etc. The gauge links in Eq. (1) are taken in the adjoint and fundamental representations, respectively,

[TABLE]

where $A_{+}=A_{+}^{a}T^{a}$ and in the adjoint representation $T^{a}_{bb^{\prime}}=if^{bab^{\prime}}$. In Eq. (2) we have introduced another notation that will be used throughout this work:

[TABLE]

The LROs (1) can be viewed as the generating functions for local quark and gluon operators

[TABLE]

where $\chi=q,g$ and one can choose, for example,

[TABLE]

where $C_{n}^{(\lambda)}(x)$ are Gegenbauer polynomials. For the $C$-parity-even operators considered here the sum goes, obviously, over odd $n$. The “coefficient functions” $\Psi^{(\chi)}_{nk}(z_{1},z_{2})$ are homogeneous polynomials of two variables $z_{1},z_{2}$ of degree $k$ for quarks and $k-1$ for gluons. The factor 6 in the definition of the gluon operator is inserted for convenience. It ensures uniform normalization of the coefficient functions later on, see Eq. (65).

Renormalized LROs are defined by the same expression with bare local operators replaced by the renormalized ones (we always assume dimensional regularization and minimal subtraction):

[TABLE]

Let us stress that the polynomials $\Psi_{nk}$ in Eq. (6) are exactly the same as in Eq. (4). Note also that the LRO on the l.h.s. of Eq. (6) does not depend on the particular choice of the local operator basis. Going over to a different basis would yield different expressions for the polynomials $\Psi_{nk}$ such that the sum is unaffected.

The renormalized local operators (5) satisfy the renormalization group equation (RGE) which has the usual form

[TABLE]

where $M$ is the renormalization scale

[TABLE]

with

[TABLE]

The anomalous dimension matrix $\gamma^{\chi\chi^{\prime}}_{nn^{\prime}}$ has a triangular form in $n$-space: its elements are nonzero only for $n\geq n^{\prime}$. The diagonal entries $\gamma_{nn}^{\chi\chi^{\prime}}$ are $2\times 2$ matrices in the $\chi,\chi^{\prime}\in\{q,\,g\}$ space, they are known to three-loop accuracy Vogt:2004mw . The off-diagonal entries in $n$-space vanish to one-loop accuracy for the special choice of local operators in Eq. (5). Beyond one loop the non-diagonal entries are non-zero and their calculation to two-loop accuracy is the main topic of this study.

The RGEs for local operators of different dimension can be combined to the RGE for their generating function, the LRO. In this representation one obtains an integro-differential equation

[TABLE]

where $\mathbb{H}(a)$ (evolution kernel) is an integral operator which has a perturbative expansion

[TABLE]

The one-loop kernel $\mathbb{H}^{(1)}_{\chi\chi^{\prime}}$ takes the form Balitsky:1987bk

[TABLE]

Here we introduced the shorthand notations for the integrals

[TABLE]

The one-loop kernel commutes with the canonical generators of the (collinear) conformal transformations,

[TABLE]

where

[TABLE]

with

[TABLE]

The expressions for the one-loop kernel in Eq. (2) have been obtained in Ref. Balitsky:1987bk by direct calculation. But in fact, they are determined completely by the anomalous dimensions $\gamma^{(1),\chi\chi^{\prime}}_{n}$, see Ref. Bukhvostov:1985rn and a discussion in Ref. Braun:2017cih . Starting from two loops this property is lost: The evolution kernel does not commute any longer with the canonical generators and cannot be restored from the anomalous dimensions alone. Nevertheless, one can simplify the calculation considerably by observing that QCD in noninteger $d-2\epsilon$ dimensions possesses a nontrivial critical point such that the $\beta$-function vanishes for the special value of the coupling, $\beta(a^{\ast})=0$, and the theory enjoys full conformal invariance. This property allows one to argue that Eq. (14) holds true for the full kernels in arbitrary order of perturbation theory with the appropriately modified (“deformed”) symmetry generators $S_{\alpha}$. A general technique for the calculation of quantum corrections to the generators of conformal transformations was developed in  Braun:2013tva ; Braun:2014vba ; Braun:2016qlg . The present case (flavor-singlet) is more complicated as compared to the discussion of flavor-nonsinglet operators in Ref. Braun:2016qlg because of a nontrivial mixing of gauge-invariant LROs with BRST Becchi:1975nq and Equation-of-Motion (EOM) operators. We consider this problem in more detail in the next Section.

3 Light-ray operators beyond one loop

It is well known that gauge-invariant (local) operators mix under renormalization with the BRST variations and EOM operators which are not gauge-invariant. A renormalized operator $[\mathcal{O}_{\chi,nk}]$ can be decomposed as Joglekar:1975nu ; Joglekar:1976eb ; Joglekar:1976pe ; collins_1984

[TABLE]

where the three terms on the r.h.s. are the gauge-invariant, BRST and EOM operators, respectively. The last two terms do not contribute to the correlators of gauge invariant operators and in many cases can be omitted. The (renormalized) LROs can be decomposed in a similar manner and one can expect that in most cases only the gauge-invariant contributions will prove to be relevant,

[TABLE]

By definition, the variation of the LRO under collinear conformal transformations $\delta_{\alpha}$, $\alpha=0,\pm$ (or, equivalently, dilatation and rotation in the $(n,\bar{n})$ plane, $D-M_{n\bar{n}}$, translation along the light-cone $P_{+}$ and special conformal transformation $K_{-}$, cf. Appendix A) is determined by the transformation properties of local operators Braun:2013tva ; Braun:2014vba ; Braun:2016qlg

[TABLE]

where $\delta_{\alpha}\widehat{\mathcal{O}}_{\chi,nk}(0)$ can be expanded over the basis of gauge-invariant local operators with finite coefficients, 222Our notation uses the following short-hands here: The last subscript $k+\alpha$ is meant to be $k-1$, $k$ and $k+1$ for $\alpha=-,0,+$, respectively.

[TABLE]

where the ellipses stand for contributions of gauge non-invariant operators. For dilatations, the structure of such terms repeats Eq. (17) – one obtains local operators that can be presented as a BRST variation, and EOM operators. For the special conformal transformations the structure of such gauge non-invariant contributions is more complicated, and they are in general non-local. A detailed discussion and explicit expressions can be found in Ref. Braun:2018mxm where it is proven that such extra terms do not contribute to correlation functions with gauge-invariant operators. For the symmetry transformations from the Poincaré group there are no such contributions, i.e., they involve only the sum over gauge-invariant operators.

Note that the expressions for the symmetry transformations depend on the choice of the local operators but this dependence is compensated by the modification of the coefficient functions. As the result, the symmetry generators in the LRO representation are defined unambiguously and do not depend on the choice of the operator basis.

Omitting gauge non-invariant contributions and substituting Eq. (20) in Eq. (19) one can represent the result as an action on the LRO of a certain linear operator

[TABLE]

The generators $S_{\alpha,\chi\chi^{\prime}}$, $\alpha=\pm,0$ are integro-differential operators in $z_{1},z_{2}$, which are $2\times 2$ matrices in the quark-gluon space, $\chi,\chi^{\prime}\in\{q,\,g\}$. From general considerations Braun:2013tva ; Braun:2016qlg it follows that they can be written in the form

[TABLE]

where the canonical generators, $S^{(0)}_{\alpha}$, are defined in Eq. (2). Thus quantum corrections to the classical symmetry generators involve the evolution kernel $\mathbb{H}$ and the operator $\Delta$ which we will refer to as conformal anomaly. Both of them are matrices in the quark-gluon space.

3.1 Flavor-singlet conformal anomaly $\Delta$

The leading-order (one-loop) conformal anomaly for the case under consideration

[TABLE]

has been calculated in Ref. Belitsky:1998gc from the analysis of the scale and conformal Ward identities. Our approach is in principle similar but differs from  Belitsky:1998gc in details. The calculation is described in Appendix B. We obtain

[TABLE]

This expression is in agreement with Ref. Belitsky:1998gc .

The general procedure how one can restore the evolution kernels from the spectrum of anomalous dimensions at the same order in perturbation theory and the conformal anomaly at one order less is explained in detail in  Braun:2017cih and it can be used for the flavor-singlet operators as well.

The central observation is that the evolution kernel at critical coupling must commute with the generators of conformal transformations,

[TABLE]

where both operators are now $2\times 2$-matrices. Expanding this equation in powers of the coupling yields a set of nested equations, e.g. up to $\mathcal{O}(a_{\ast}^{2})$,

[TABLE]

Eq. (26a) is a first-order homogeneous differential equation and its general solution (which we call “invariant operator” in what follows) takes the generic form

[TABLE]

where $h_{\rm inv}(\tau)$ are functions of $\tau=\frac{\alpha\beta}{\bar{\alpha}\bar{\beta}}$ dubbed “the conformal ratio”. Indeed, the one-loop evolution kernel Eq. (2) can be rewritten in this form Braun:2014vba . Eq. (26b), in turn, can be viewed as an inhomogeneous differential equation for the two-loop evolution kernel $\mathbb{H}^{(2)}$ and is solved by Braun:2017cih

[TABLE]

where $\mathbb{H}_{\rm inv}^{(2)}$ can be any invariant operator that takes the form Eq. (27) and the operators $\mathbb{X}^{(1)}$ and $\mathbb{T}^{(1)}$ are $2\times 2$-matrices which must satisfy the constraints

[TABLE]

Solving these two equations we get

[TABLE]

The operator $\mathbb{T}^{(1)}$ takes the form

[TABLE]

This expression differs from the evolution kernel in Eq. (2) only by the insertion of a factor $\ln(1-\alpha-\beta)$. Note that the constraints in Eq. (29) fix $\mathbb{X}^{(1)}$ and $\mathbb{T}^{(1)}$ only up to contributions of invariant operators. Such possible extra terms enter the two-loop kernel $\mathbb{H}^{(2)}$ through the product with $\mathbb{H}^{(1)}$ which is itself an invariant operator. One can easily verify that the product of two invariant operators Eq. (27) is again an invariant operator. Therefore, the freedom in the definition of $\mathbb{X}^{(1)}$ and $\mathbb{T}^{(1)}$ corresponds to a redefinition of the two-loop invariant evolution kernel $\mathbb{H}_{\rm inv}^{(2)}$ which still has to be determined.

The expressions given above define the non-invariant part of the two-loop evolution kernels in the factorized form

[TABLE]

as a product of relatively simple integral operators (3.1), (3.1). For certain applications and in particular for numerical studies of the scale dependence of the GPDs, it can be advantageous to have explicit expressions for these products.

The results can be written in the form

[TABLE]

where

[TABLE]

For the pure singlet-contribution to the quark kernel we get $\Delta r_{qq}^{\rm PS}=\vartheta_{qq}^{\rm PS}=0$ and

[TABLE]

and for the non-singlet contribution (cf. Ref. Braun:2014vba )

[TABLE]

For the off-diagonal kernels we obtain $\Delta r_{qg}=\vartheta_{qg}=0$,

[TABLE]

and

[TABLE]

where the “plus” distribution is defined as

[TABLE]

Finally, for the gluon-gluon kernel we obtain

[TABLE]

The expressions in Eqs. (3.1)-(3.1) supplemented by the invariant kernels which are calculated in the next Section provide one with the complete two-loop flavor-singlet evolution kernels in the LRO representation. These are our main results.

4 Anomalous dimensions vs invariant kernels

Substituting the LRO in Eq. (10) by its expansion in terms of local operators and comparing the resulting expression with the RGE (7) one gets

[TABLE]

Conformal operators of the lowest dimension for the given spin correspond to the highest weights of the representation and are annihilated by $S_{-}=S_{-}^{(0)}$. As a consequence, the coefficient functions of the operators with $k=n$ are translation-invariant and by dimension counting $\Psi_{nn}^{q}\sim z_{12}^{n}$ and $\Psi^{g}_{nn}\sim z_{12}^{n-1}$. It follows that

[TABLE]

where the coefficients $h_{\chi\chi^{\prime}}$ are related to the matrix of anomalous dimensions $\gamma_{\chi\chi^{\prime}}(n)$ in the normalization chosen in Ref. Vogt:2004mw as

[TABLE]

Eq. (26b) allows one to find the spectrum of the invariant kernel $\mathbb{H}_{inv}^{(2)}$ from the known results for the two-loop anomalous dimensions. Following the notation in Ref. Vogt:2004mw we split the two-loop flavor-singlet quark kernel in two parts – flavor non-singlet and pure-singlet, $\mathbb{H}^{(2)}_{qq}=\mathbb{H}^{(2)}_{\rm NS}+\mathbb{H}^{(2)}_{\rm PS}$. The expression for the non-singlet part can be found in Ref. Braun:2017cih , while for the pure singlet one finds

[TABLE]

where $h^{(2)}_{\rm PS}(n)$ and $h^{(1)}_{\chi\chi^{\prime}}(n)$ are related to the corresponding anomalous dimensions as in Eq. (43), $\partial_{n}=\frac{\partial}{\partial n}$ and

[TABLE]

Using

[TABLE]

and the two-loop expression for the pure-singlet anomalous dimension from Ref. Vogt:2004mw , we get

[TABLE]

The corresponding kernel as a function of the conformal ratio

[TABLE]

takes the form

[TABLE]

In a similar manner one obtains from Eq. (26b) the equations for the remaining entries:

[TABLE]

The relevant eigenvalues of the $\mathbb{X}$ kernels are equal to

[TABLE]

cf. Eq. (45) for $X^{(1)}_{gq}$, and the eigenvalues for the diagonal entries of the one-loop evolution kernels take the familiar form

[TABLE]

Their derivatives can be written as

[TABLE]

Combining all these expressions we obtain, after some algebra, the following results for the eigenvalues of the invariant kernels

[TABLE]

[TABLE]

and

[TABLE]

It remains to restore the invariant kernels in the LRO representation from these results for the spectrum. The general expression for $\mathbb{H}^{\text{inv}}$ can be written as a $2\times 2$ matrix

[TABLE]

The entries are invariant operators which we parameterize as follows

[TABLE]

where $\tau=\alpha\beta/\bar{\alpha}\bar{\beta}$. Adding the pure-singlet kernel (49) to the expression for the flavor non-singlet kernel given in Ref. Braun:2014vba ; Braun:2017cih one obtains

[TABLE]

where $\frac{1}{N_{c}}=C_{A}-2C_{F}$ for an SU$(N_{c})$ gauge theory. For the gluon-quark kernel one obtains using the expressions collected in Appendix D,

[TABLE]

For the quark-gluon kernel we find $\Gamma_{qg}=r_{qg}=0$ and

[TABLE]

Finally, for the gluon-gluon kernel we obtain

[TABLE]

and

[TABLE]

We close this discussion with a remark on the so-called reciprocity symmetry of the invariant kernels Basso:2006nk . This symmetry arises, technically, from the observation Basso:2006nk that invariant kernels can in general be presented in terms of the quadratic Casimir operator of the collinear conformal group. For an operator with conformal spin $j$ (in our case $j=n+2$) and its anomalous dimension given by $\gamma(j)=f(j+\frac{1}{2}\gamma(j))$, the asymptotic expansion of the function $f(j)$ for large $j$ consists of terms invariant under $j\to 1-j$. It was shown in Ref. Braun:2017cih that the function $f(j)$ defined by this equation gives the eigenvalues of the invariant kernel, $\mathbb{H}^{\text{inv}}$, hence the eigenvalues of $\mathbb{H}^{\text{inv}}$ should have (and indeed they have) the corresponding invariance property. The reciprocity relation has been checked on many examples for the situations where only one operator exists for a given conformal spin. In the situation that there are two and or more operators of the same conformal spin, the reciprocity cannot be expected in general for the off-diagonal elements of the anomalous dimension matrix, because they depend explicitly on the assumed normalization for the operators Basso:2006nk . Since the evolution kernels for the LROs do not depend on the basis of local operators, one should expect, however, that in this representation the reciprocity holds for the off-diagonal elements (kernels) as well. Indeed, one can verify that the eigenvalues of the off-diagonal invariant kernels, Eqs. (4) and (4) are invariant under $j\to 1-j$ with $j=n+2$. For the invariant kernels themselves this condition implies an expansion $\sum_{k}h_{k}(\tau)\ln^{k}\tau$ for $\tau\to 0$, where $h_{k}(\tau)$ is an analytic function in the vicinity of $\tau=0$.

5 Anomalous dimension matrix for local operators

For a comparison with Ref. Belitsky:1998gc and also for the application to the scale-dependence of the meson light-cone distribution amplitudes, it is desirable to have the results also in a different form, as an anomalous dimension matrix for local operators in the Gegenbauer basis (5).

The coefficient functions $\Psi^{(\chi)}_{nk}(z_{1},z_{2})$ in the expansion of the quark and gluon LROs over this basis (5)

[TABLE]

can be obtained by the repeated application of the canonical $S^{(0)}_{+}$ generator to the coefficient function of the highest-weight operator $k=n$,

[TABLE]

where

[TABLE]

These polynomials are mutually orthogonal and form a complete set of functions w.r.t. the canonical $SL(2)$ scalar product (see, e.g., Ref. Braun:2011dg )

[TABLE]

where

[TABLE]

The local operators (5) can be obtained by the projection of the LROs on the corresponding coefficient function

[TABLE]

In order to translate the action of the evolution kernel $\mathbb{H}$ in the LRO representation in Eq. (28) in terms of local operators we adopt the following dictionary:

[TABLE]

where $\mathbb{A}\in\{\mathbb{H},\mathbb{X},\mathbb{T}\}$. Note that for all operators from this set $[\mathbb{A},S_{-}]=0$ and, therefore, we can choose $k=n$ in Eq. (5) for convenience without loss of generality. Moreover, it is worth to be mentioned that for any invariant operator $\mathbb{A}_{\rm inv}^{\chi\chi^{\prime}}$ the matrix $\bm{A}_{\rm inv}^{\chi\chi^{\prime}}$ is diagonal:

[TABLE]

Making use of these properties we find the following results for the local matrix elements of the $\mathbb{X}$ operators (3.1)

[TABLE]

with the eigenvalues $X^{(1)}_{\chi\chi^{\prime}}(n)$ given in Eqs. (45) and (51), the discrete step function defined as

[TABLE]

and the following matrices:

[TABLE]

where

[TABLE]

For the operator $\mathbb{T}^{(1)}$ defined in Eq. (3.1) we obtain

[TABLE]

with the one-loop anomalous dimensions in Eqs. (46), (4), (4).

This yields the following representation for the anomalous dimension matrix $\bm{H}_{nm}^{(2)\chi\chi^{\prime}}$ of the local operators

[TABLE]

The eigenvalues $h_{\chi\chi^{\prime}}^{\rm inv}(n)$ of the two-loop invariant kernels can be found in Eqs. (47), (4), (4) and (4).

The result in this form can directly be compared to the two-loop anomalous dimension matrix calculated in Ref. Belitsky:1998gc . The diagonal elements $\bm{H}^{(2)\chi\chi^{\prime}}_{nn}$ for each channel reproduce the well-known eigenvalues $h_{\chi\chi^{\prime}}^{(2)}(n)$ of the evolution kernel by construction, Eqs. (44), (4), and in both calculations coincide (up to a different normalization) with the anomalous dimensions from Ref. Vogt:2004mw . For the non-diagonal elements note that our definition of the gluonic local operator in Eq. (5) differs from the one used in Ref. Belitsky:1998gc by a factor 6, i.e.

[TABLE]

and also that the anomalous dimensions in Ref. Belitsky:1998gc are written as an expansion in $\alpha_{s}/(2\pi)$. Taking these differences into account we find a perfect agreement.

6 Evolution kernels in the momentum fraction representation

Evolution equations in the momentum fraction representation can be derived straightforwardly from the kernels in the coordinate representation. Let $f(x_{1},x_{2})$ be a Fourier transform of the position-space distribution $f(z_{1},z_{2})$ defined as

[TABLE]

The general expression for an evolution kernel in the momentum fraction representation takes the form

[TABLE]

where $\mathcal{D}u=du_{1}du_{2}\delta(x_{1}+x_{2}-u_{1}-u_{2})$ and the $\delta$-function is due to the momentum conservation. Following the decomposition in Eq. (4) we split the kernel $\mathcal{H}$ into three parts

[TABLE]

where $\widehat{\mathcal{H}}_{\vartheta}$, $\mathcal{H}_{\omega}$ originate from the contributions which have the form in the coordinate space (LRO) representation

[TABLE]

It is easy to show that the corresponding kernels in the momentum fraction representation are given by the following expressions:

[TABLE]

The $\Theta$ function is defined as follows

[TABLE]

and functions $A_{\omega},B_{\omega},C_{\omega}$ are given by the following integrals

[TABLE]

where it is implied that $x_{1}+x_{2}=u_{1}+u_{2}$ and

[TABLE]

The kernel $\omega(\alpha,\beta)$ is a single-valued function in the simplex $0<\alpha+\beta<1$. It is easy to see that if the variables $x_{i},u_{i}$ belong to the regions determined by the corresponding $\Theta$-functions, then the arguments $\alpha_{x},\beta_{x}$ lie inside the simplex and the corresponding integrals are unambiguously defined. In the present case, $\omega(\alpha,\beta)$ is symmetric under the interchange of its two arguments $\omega(\alpha,\beta)=\omega(\beta,\alpha)$ (this is not always the case, see Ref. Braun:2009mi ), and as a consequence $C_{\omega}(x_{1},x_{2},u_{1},u_{2})=B(x_{2},x_{1},u_{2},u_{1})$. Provided that the functions $B_{\omega},C_{\omega}$ are continued in an appropriate way beyond their analyticity domain one gets $A_{\omega}=B_{\omega}-C_{\omega}$. Note also that all these kernels are effectively functions of two variables, $x=x_{1}/(x_{1}+x_{2})$ and $y=u_{1}/(x_{1}+x_{2})$, e.g. $A_{\omega}=a(x,y)/(x_{1}+x_{2})$, etc.

All integrals that arise in the transition from the position to the momentum fraction representation can be taken analytically in terms of elementary functions and $\text{Li}_{2}$ polylogarithms. The resulting expressions are collected in Ref. Belitsky:1999hf and are implemented in the NLO evolution FORTRAN code by Freund and McDermott Freund:2001rk ; Freund:2001hd  333http://durpdg.dur.ac.uk/hepdata/dvcs.html. We have checked numerically that our results for the kernels in the momentum fraction representations agree with the kernels implemented in this code for all color structures. In this way we can confirm the results of Ref. Belitsky:1999hf , where these kernels are given in analytic form 444See footnote 17 in Ref. Freund:2001hd ..

7 Summary

We have presented a re-derivation of the two-loop flavor-singlet evolution equations for the leading-twist operators in off-forward kinematics, based on using conformal symmetry of QCD in non-integer $d-2\epsilon$ space-time dimensions at the critical point. This case is more complicated as compared to the flavor-nonsinglet evolution, both, technically and conceptually, due to potentially dangerous contributions of gauge-noninvariant operators. Our analysis is based on studies of conformal Ward identities in Ref. Braun:2018mxm , where it is proven that such extra terms do not contribute to correlation functions with gauge-invariant operators.

The results in this work are given in the position-space or light-ray operator representation. This form has some technical advantages and can also be interesting in applications to lattice QCD calculations. Expanding our results in powers of the field separation we reproduce the results for the mixing matrices for flavor-singlet local operators derived in Ref. Belitsky:1998gc . The evolution kernels in momentum fraction space can be obtained from our expressions by simple integration. We have verified numerically that the resulting kernels agree with the kernels that are implemented in the NLO GPD evolution code by Freund and McDermott Freund:2001rk ; Freund:2001hd which is based on the analytic expressions derived in Ref. Belitsky:1999hf .

To summarize, the two-loop evolution equations for generalized parton distributions have now been derived independently by two groups, and perfect agreement is found. In view of the projected high statistics of the relevant experiments at the JLAB 12 GeV upgrade and, in future, the EIC, it is necessary to implement and use these results now in NLO analyses of deeply-virtual Compton scattering and similar reactions.

Acknowledgments

We are grateful to V. Guzey for providing us with the source code for the NLO GPD evolution by Freund and McDermott Freund:2001rk ; Freund:2001hd . This work was supported by the DFG with the grants MO 1801/1-3 (A.M.) and SFB/TRR 55 (M.S.).

Appendices

Appendix A Scale and conformal transformations

Scale ($D$) and conformal ($K$) field transformations for the fundamental fields take the generic form

[TABLE]

in particular

[TABLE]

where $\Delta_{\Phi}=\dim\Phi$ are the canonical dimensions of the fields. It is convenient to choose them in $4-2\epsilon$ dimensions to be exactly the same as in the four-dimensional theory,

[TABLE]

For this choice the gluon strength tensor $F_{\sigma\rho}$ transforms in a covariant way

[TABLE]

and the covariant derivative of the ghost field $D_{\nu}c$ transform as a vector field,

[TABLE]

The variation of the different parts of the QCD action under the special conformal transformation takes the form

[TABLE]

Note that the ghost and the gauge-fixing terms break conformal symmetry explicitly even in $d=4$ dimensions. Summing up all contributions one obtains

[TABLE]

Here $\mathcal{J}_{\rho}(x)=\bar{q}(x)\gamma_{\rho}q(x)$ is the conserved flavor-singlet current and

[TABLE]

is a BRST variation, $\mathcal{B}_{\mu}=\delta_{\rm BRST}(\bar{c}A_{\mu})$.

Appendix B One-loop conformal anomaly

Our analysis follows closely the lines of Ref. Braun:2016qlg so that in this Appendix we concentrate mainly on specific problems that arise for the flavor-singlet operators.

Consider Braun:2016qlg the correlation function of two (renormalized) LROs stretched in different directions

[TABLE]

where $z=\{z_{1},z_{2}\}$, $w=\{w_{1},w_{2}\}$ and $n$, $\bar{n}$ are two auxiliary light-like vectors, $n^{2}=\bar{n}^{2}=0$. In what follows we assume that $(n\bar{n})=1$ and $(x\cdot n)=(x\cdot\bar{n})=0$. Starting from the representation in terms of local conformal operators it can be shown, see Ref. Braun:2016qlg , that the correlation function (B.96) satisfies the constraint

[TABLE]

where the operator $S_{+,\sigma\chi}$ is the generator of special conformal transformations defined in Eq. (21). The explicit expression for $S_{+,\sigma\chi}$ can be found from the analysis of the conformal Ward identity (CWI) for the correlation function (B.96). Making the corresponding change of variables in the path integral representation for this correlation function one obtains

[TABLE]

where $\delta$ stands for the conformal variation along the $\bar{n}$ direction, $\delta=\bar{n}^{\mu}\delta_{K^{\mu}}$. The second term on the l.h.s. is solely responsible for the r.h.s. of Eq. (B.97) while the two others contribute to the l.h.s. of this relation. The first of them, the variation of the LRO operator, takes the form

[TABLE]

The renormalization factor $Z$ in this equation is an integral operator in $z_{1},z_{2}$ and it does not commute with $S_{\chi,+}^{(0)}$. However, since $S_{\chi,+}^{(0)}$ does not have any $\epsilon$ dependence,

[TABLE]

Since all singular terms must cancel in the final result we can drop them.

The second term that contributes to the l.h.s. of Eq. (B.97) is due to the conformal variation of action as given in Eq. (A.94). As the first step, it has to be re-expanded in terms of renormalized operators. In Landau gauge that we will assume from now on, the corresponding expression takes the form

[TABLE]

where $\mathcal{B}^{\mu}$ is defined in Eq. (A.95), $\mathcal{L}^{\prime}_{R}(x)=\mathcal{L}_{R}(x)-\frac{1}{2}\,\partial^{\rho}[\mathcal{J}_{\rho}](x)$, $\Omega_{\Phi}=\Phi{\delta S_{R}}/{\delta\Phi}$, $\Omega_{q\bar{q}}=\Omega_{q}+\Omega_{\bar{q}}$. Note that the anti-ghost EOM term does not contribute to the correlation function in question since the LROs $\widehat{\mathcal{O}}_{\chi}$ do not contain (anti)ghost fields. The terms $\Omega_{A}$ and $\partial_{\mu}[\mathcal{B}^{\mu}]$ contain contributions of the gauge parameter proportional to $1/\xi$, which, however, cancel each other. For nonzero $\xi$ the correlation function of $\mathcal{B}_{\mu}$ with any gauge invariant functional $\mathcal{F}$ vanishes, implying that

[TABLE]

Integrating by parts the contributions of the EOM terms in $\langle{\delta S_{R}\,\widehat{\mathcal{O}}_{\chi}^{(n)}(0,z)\,\widehat{\mathcal{O}}_{\chi^{\prime}}^{(\bar{n})}(x,w)}\rangle$ we can write, e.g., for the first operator insertion,

[TABLE]

where we introduce the notation $\mathcal{D}_{\Phi}(x)=\Phi(x)\frac{\delta}{\delta\Phi(x)}$, $\mathcal{D}_{q\bar{q}}(x)=\mathcal{D}_{q}(x)+\mathcal{D}_{\bar{q}}(x)$.

Let us consider first the quark operator, $\widehat{\mathcal{O}}_{q}(0,z)$. In order to streamline the notation, hereafter, we do not display the dependence on the auxiliary vector $n$ and the space-time coordinates, if this is clear from the context. Our goal is to calculate the one-loop correction to the generators. It can be shown that there are no BRST or EOM counterterms to the quark LRO at one loop, i.e. $\widehat{\mathcal{O}}_{q}(z)=Z_{qq}\mathcal{O}_{q}(z)+Z_{qg}\mathcal{O}_{g}(z)+O(a^{2})$. It is also obvious that the one-loop correction to the quark part of the generator, $S_{qq}$, is the same as in the non-singlet case. We consider therefore only the off-diagonal entry $S_{qg}$, i.e. we are looking for the contributions of the gluon LRO $\widehat{\mathcal{O}}_{g}$ on the r.h.s. of Eq. (B.103).

First, note that the term $2\gamma_{c}\mathcal{D}_{A}(x)\widehat{\mathcal{O}}_{q}(z)$ gives rise to a $O(a^{2})$ contribution only and can be omitted. Second, there are no gluon pair counterterms for the product $\left[\mathcal{L}^{YM}(x)\right]\widehat{\mathcal{O}}_{q}(z)$ at one loop. Thus, the only relevant contribution comes from the quark EOM term. We write

[TABLE]

The gluon contribution of interest corresponds to the last term in this expression. Taking into account that the renormalization factor is related to the evolution kernel as $Z^{-1}_{qg}=-\frac{a}{2\epsilon}\mathbb{H}_{qg}^{(1)}+O(a^{2})$ and multiplying this relation by $\epsilon-\gamma_{q}$, one can bring this contribution to the form

[TABLE]

This result implies that $\Delta^{(1)}_{qg}=0$ in agreement with Ref. Belitsky:1998gc .

Next, consider the $S_{gq}$ generator. To this end we again omit the $\mathcal{D}_{A}$ term in Eq. (B.103) since it is $O(a^{2})$. The quark EOM contribution in Eq. (B.103) can be handled as above yielding

[TABLE]

In addition we have to consider the contribution from the product $\left[\mathcal{L}^{YM}(x)\right]\widehat{\mathcal{O}}_{g}(z)$. This product is divergent when $x\to 0$ so that one has to add and subtract the corresponding pair counterterms. The counterterm $\sim\mathcal{O}_{q}$ corresponds to the Feynman diagram shown in Fig. 1

where the crossed circle stands for the insertion of the vertex

[TABLE]

The terms involving $(\partial A)$ vanish since the propagators are transverse in Landau gauge. The gluon line with the crossed-circle insertion (modified propagator) can be written in the form

[TABLE]

where $G_{\mu\nu}$ is the usual gluon propagator in Landau gauge and

[TABLE]

Taking into account the contribution of the crossing-symmetric diagram one obtains (divergent part only) 555The contributions due to $H_{\mu\nu}$ vanish in the sum of all diagrams.

[TABLE]

where $\{\ldots,\ldots\}_{+}$ stands for the anti-commutator. Multiplying this expression by $(-\beta(a)/2a)$ and adding Eq. (B.106) we obtain the result for $\Delta^{(1)}_{gq}$ given in Eq. (24).

The calculation of the pure gluon contribution, $\Delta_{gg}$, is more involved. It was shown in Braun:2018mxm that the r.h.s. of Eq. (B.103) can be represented as sum of the fully renormalized product $[F_{\mu\nu}^{2}(x)\mathcal{O}_{g}(z)]$ which enters with the coefficient $\beta(a)/a$, contributions of gauge-invariant operators and gauge-noninvariant contributions that are BRST variations of nonlocal operators or nonlocal EOM operators which drop out of all correlation functions with gauge-invariant operators, hence also of Eq. (B.98). Nevertheless, we need to know the structure of the gauge-noninvariant terms in order to correctly separate the gauge-invariant contributions. Consider first the BRST and EOM counterterms to the gluon LRO $\mathcal{O}_{g}$. Since they should be twist-two operators, there are two possible structures:

[TABLE]

Here we tacitly assume that the color indices are contracted in some way. Both operators contain singular terms in $\xi$ that have to cancel in the sum. Therefore they can appear only in the combination, $B_{k}(z_{1},\ldots,z_{k})-E_{k}(z_{1},\ldots,z_{k})$. The lowest, $k=2$, counterterm can be restored by calculating the one-loop $F_{+\mu}F_{+\mu}\mapsto\bar{c}c$ diagram. The result reads

[TABLE]

The BRST variation of this operator is given by a sum of EOM contributions which, in agreement with the general statement, do not contain the anti-ghost EOM. The contributions for $k=3,4,\ldots$ can be neglected to our accuracy as they are higher order in the coupling.

Next we notice that for the calculation of $\Delta_{gg}$ to one-loop accuracy one can omit the quark EOM term $\sim\mathcal{D}_{q\bar{q}}$ in Eq. (B.103), and the remaining contribution $\sim\mathcal{D}_{A}(x)$ can be written as

[TABLE]

Here

[TABLE]

where the gauge links (in the adjoint representation) are tacitly implied. Since Eq. (B.113) is multiplied by the ghost anomalous dimension $2\gamma_{c}$ the contribution of $\mathcal{O}_{3g}$ is higher order in the coupling and can be omitted to the $O(a)$ accuracy. Thus the first term in Eq. (B.103) can be presented in the form

[TABLE]

where $(B-E)_{2}$ is given in Eq. (B.112) and the ellipses stand for $k>2$ (gauge-noninvariant) counterterms which contribute only starting from $O(a^{2})$. The term with the BRST variation $B_{2}$ also drops out from the correlation function at order $O(a)$. The remaining EOM term can be rewritten (inside the correlation function) with the required accuracy as

[TABLE]

This is the only contribution from the BRST and EOM counterterms at this order.

Note that the operator in Eq. (B.117) is not a BRST variation so that, according to the general statement of Ref. Braun:2018mxm , this contribution must be complemented by contributions of other diagrams to produce a BRST variation or EOM operator in the sum, which is not easy to see on a diagrammatic level. For our purposes it is, however, sufficient to note that all such “unwanted” contributions contain either $A_{+}$ or $A_{-}$ fields, and can be avoided if one looks for the contributions of transverse gluons.

In order to calculate the pair counterterm to the operator product $\left[\mathcal{L}^{YM}(x)\right][{\mathcal{O}}_{g}(z)]$ it is convenient to make some rearrangements. Namely, we rewrite

[TABLE]

The integral of the first piece can be written as (up to $(\partial A)$ contributions)

[TABLE]

where $\mathcal{L}_{3,4}^{YM}$ stand for the three- and four-gluon vertices. The ghost terms can be omitted as they do not contribute to the gauge-invariant counterterms. In the second piece, the BRST operator does not contribute to the correlation function and the $\Omega_{A}$-term gives

[TABLE]

All these contributions are multiplied by ${(-\beta(a)/2a)}=2(\epsilon-\gamma_{g})$.

The main contribution to the CWI is due to the pair counterterms to

[TABLE]

The corresponding diagrams are shown in Fig. 2 (adding the symmetric contributions is implied).

They can be related to the counterterm diagrams for the gluon LRO with some extra terms.

Since the four-gluon vertex does not contain derivatives, the counterterm corresponding to the last diagram in Fig. 2 can be rewritten as

[TABLE]

where $G_{\mu\nu}^{ab}(y,z_{1},z_{2})$ is the usual QCD diagram with the four-gluon vertex and $\overline{G}_{\mu\nu}^{ab}(y,z_{1},z_{2})$ is the diagram with the modified vertex, $\overline{G}_{\mu\nu}^{ab}(y,z_{1},z_{2})\equiv 2(\bar{n}y)G_{\mu\nu}^{ab}(y,z_{1},z_{2})$.

For the rest of the contributions in Fig. 2 that involve a three-gluon vertex (which contains a derivative), this relation has to be modified as follows. Let $y_{1}$ be the coordinate of the modified vertex (insertion of $(\bar{n}\cdot y_{1})$). Then

[TABLE]

where we suppressed color and Lorentz indices on the gauge field. In this expression, as above, $G(y_{1},y_{2},z_{1},z_{2})$ is the diagram with the standard QCD three-gluon vertex, $\overline{G}(y_{1},y_{2},z_{1},z_{2})=\ (\bar{n}y_{1})G(y_{1},y_{2},z_{1},z_{2})$ and $\widetilde{G}_{1}(y_{1},y_{2},,z_{1},z_{2})$ is the Feynman diagram obtained from $G(y_{1},y_{2},z_{1},z_{2})$ by the replacement of the triple vertex $V(q,r,k)$ at the position $y_{1}$ by the new vertex $\widetilde{V}(q,r,k)$ defined as

[TABLE]

Using these relations the divergent parts of the contributions in Fig. 2 can be presented in the form

[TABLE]

Here $\text{CT}(O_{g})$ stands for the counterterms to the operator $\mathcal{O}_{g}$, $\widetilde{D}_{i}$ are the first four diagrams in Fig. 2 (with the triple-gluon vertex) with the crossed vertex replaced by Eq. (B.124) and $D^{\prime}_{i}$ are the diagrams shown in Fig. 3 (upper row), where the big black circle stands for the vertices generated by the operator $(n\bar{n})A_{\mu}(z_{2})-\bar{n}^{\mu}A_{+}(z_{2})+(n\bar{n})z_{2}F_{+\mu}(z_{2})\,$.

The first contribution in Eq. (B.125) absorbs all terms with the derivative $\mathcal{D}_{A}$ in Eqs. (B) and (B.122) but also generates additional contributions, the diagrams $D^{\prime}_{i}$ shown in the upper row in Fig. 3, that have, therefore, to be subtracted.

We recall that Eq. (B.125) has to be multiplied by factor ${(-\beta(a)/2a)}=\epsilon-\gamma_{g}$ and also the contribution due to Eq. (B.120) has to be added, which reads

[TABLE]

The first term in Eq. (B.125) together with Eq. (B.126) gives

[TABLE]

where the operators $\mathcal{O}_{2,g}(z)$ and $\mathcal{O}_{3,g}(z)$ are defined in Eqs. (B.114) and (B.115), respectively. This expression has still to be rewritten as the sum of renormalized operators minus the corresponding counterterms. Since $\gamma_{g}(a_{*})=\epsilon$ the contributions of renormalized operators vanish at the critical point. For the first term we can write

[TABLE]

Combining this expression with Eq. (B.113) and taking into account that in Landau gauge $2\gamma_{c}=-\gamma_{A}-\gamma_{g}$ one obtains the usual contribution to the generator $S_{+}$

[TABLE]

The counterterms to the operator $\mathcal{O}_{3g}$ in Eq. (B.127) come from the diagrams shown in the middle row in Fig. 3. The first of them cancels with the first diagram in the upper row. The remaining two diagrams in the upper and the middle row do not cancel completely but their sum can be simplified to the diagrams shown in lowest row. Strictly speaking, there is also a contribution of the form $(z_{1}-z_{2})\mathcal{O}_{g}(z_{1},z_{2})$ that cancels in the sum with the symmetric diagram ($z_{1}\leftrightarrow z_{2}$). Finally, the counterterm diagrams for the operator $\mathcal{O}_{2g}$ are shown in Fig. 4.

The third and the fourth diagram in Fig. 4 cancel with the first and the second diagrams in the lowest row in Fig. 3, respectively. The remaining ones give rise to the conformal anomaly $\Delta_{gg}$. These are:

• •

The first four diagrams in Fig. 2 with the triple gluon vertex (B.124).

• •

The last diagram in the lowest row in Fig. 3.

• •

The first, second, fifth and sixth diagram in Fig. 4.

Technically, it is convenient to combine the first two diagrams in Fig. 2 with the first two diagrams in Fig. 4 as there are strong cancellations. The calculation of all other diagrams is rather straightforward. Note, that the final answer contains not only the term $F_{\mu+}F_{\mu+}$ which we are interested in, but also gauge non-invariant contributions. These terms contain the field $A_{+}$ or $A_{-}$ or both of them and must be part of the BRST and EOM operators which can appear on the r.h.s. of Eq. (B.103). We did not keep track of such contributions. But we have checked that the gauge-noninvariant terms, $(n\bar{n})A_{\mu}(z_{1}n)F_{+\mu}(z_{2}n)$, which are neither BRST nor EOM operators, cancel out in the sum of all diagrams. The resulting expression for the conformal anomaly is presented in Eq. (24). It has a very simply form and coincides with the result of Belitsky and Müller Belitsky:1998gc , that is the ${}^{\chi\chi^{\prime}}\mathcal{K}^{\omega}$ kernels in their notation.

Appendix C RG identities in Landau gauge

The analysis of the Ward identities in this work is based on the re-expansion of the renormalized Lagrangian $2\epsilon\mathcal{L}_{R}$ in terms of renormalized operators. The corresponding general expression in arbitrary covariant gauge reads Spiridonov:1984br ; Belitsky:1998gc ; Braun:2016qlg

[TABLE]

Here $\mathcal{L}^{\prime}_{R}=\mathcal{L}_{R}-\frac{1}{2}\partial^{\rho}(\bar{q}_{0}\gamma_{\mu}q_{0})$, $\gamma_{\Phi}$ are the anomalous dimensions of the fields, $\Omega_{\Phi}(x)=\Phi(x){\delta S_{R}}/{\delta\Phi(x)}$ and $\Omega_{\mu}=\bar{c}D_{\mu}c-\partial_{\mu}\bar{c}c$ is a conserved current, $\partial_{\mu}[\Omega^{\mu}]=\Omega_{c}-\Omega_{\bar{c}}$. The coefficients $z_{b}(g,\xi)$ and $z_{c}(g,\xi)$ are at this stage unknown functions.

The gauge fixing term, $[\mathcal{L}^{gf}]=\frac{1}{2\xi}(\partial A)^{2}$, can be written as a combination of BRST and EOM operators,

[TABLE]

An arbitrary diagram with an insertion of this operator vanishes in Landau gauge. Hence, this contribution can safely be omitted. Next, whereas all counterterms are polynomials in the gauge parameter $\xi$, the terms $\Omega_{A}$ and $\partial_{\mu}[\mathcal{B}^{\mu}]$ on the r.h.s. of Eq. (C) contain contributions $\sim 1/\xi$ which have to cancel in the sum. This implies that $z_{g}(g,\xi)=\gamma_{A}+\gamma_{g}+O(\xi)$. The last coefficient $z_{c}(g,\xi)$ can be fixed in Landau gauge by observing that in this gauge the QCD action has an additional discrete symmetry corresponding to the interchange of ghost and antighost fields, $c\to-\bar{c}$, $\bar{c}\to c$. Indeed it is easy to check that

[TABLE]

Since the last term can be omitted in Landau gauge, the action is invariant with respect to this transformation. Clearly Eq. (C) should respect this property. Taking into account that $\Omega_{c}+\Omega_{\bar{c}}$ is invariant under this transformation, $\Omega^{\mu}\mapsto-\Omega^{\mu}$, and

[TABLE]

one obtains $z_{c}(g,\xi\to 0)=\gamma_{c}=-(\gamma_{A}+\gamma_{g})/2$. Using these expressions we obtain the following result

[TABLE]

which holds in Landau gauge. Note, that here we do not keep the gauge fixing term $\mathcal{L}^{gf}$, since it does not contribute in this gauge.

Appendix D Evolution kernels vs. anomalous dimensions

In this Appendix we collect the expressions for the invariant kernels that appear as building blocks in two-loop evolution equations and the anomalous dimensions. The anomalous dimension $\gamma(j)$ corresponding to the kernel $h(\tau)$ is defined as follows:

[TABLE]

The inverse relation reads Braun:2014vba

[TABLE]

where it is assumed that all poles in anomalous dimensions are located to the left of the integration contour which is along the imaginary axis. Here $P_{j}(z)$ is the Legendre function,

[TABLE]

Below we list the relevant expressions for the anomalous dimensions $\gamma(j)$ (left column) and the corresponding ones for the kernels $h(\tau)$ (right column):

[TABLE]