Higher-twist corrections to gluon TMD factorization

TL;DR

This paper investigates the impact of higher-twist power corrections on gluon transverse momentum dependent (TMD) factorization in hadron-hadron collisions, providing insights into the accuracy of TMD descriptions.

Contribution

It introduces calculations of higher-twist power corrections to gluon TMD factorization, advancing understanding of non-leading contributions in particle production.

Findings

Quantified the size of higher-twist corrections in gluon TMD processes

Provided a framework for incorporating power corrections into TMD factorization

Enhanced the precision of gluon TMD models in high-energy collisions

Abstract

We calculate power corrections to TMD factorization for particle production by gluon-gluon fusion in hadron-hadron collisions.

Equations

$$\frac{d\sigma }{d\eta d^2{q}_{\perp }}=\sum _f\int d^2{b}_{\perp }{e}^{i(q,b{)}_{\perp }}{\mathcal{D}}_{f/A}(x_A,b_{\perp },\eta ){\mathcal{D}}_{f/B}(x_B,b_{\perp },\eta )\sigma (ff\to H)$$

$$+\mathrm{power}\mathrm{corrections}+\mathrm{Y}-\mathrm{terms}$$

$${\mathcal{L}}_{\Phi }=g_{\Phi }\int d^{4}x\Phi (x)g^2{F}^2(x),F^2(x)\equiv F_{\mu \nu }^a(x)F^{a\mu \nu }(x)$$

$$g_{H}~{}=~{}{1\over 48\pi^{2}v}\big{(}1+{11\over 4\pi}\alpha_{s}+...\big{)}$$

$$d\sigma =\frac{d^{3}q}{2E_q(2\pi {)}^3}\frac{g_{\Phi }^2}{2s}W(p_A,p_B,q)$$

$$W(p_A,p_B,q)$$

$$=$$

$$W(p_A,p_B,q)=\sum _X\int d^{4}x{e}^{-iqx}⟨p_A,p_B|g^2{F}^2(x)|X⟩⟨X|g^2{F}^2(0)|p_A,p_B⟩$$

$$=\underset{t_i\to -\infty }{\overset{t_f\to \infty }{\lim }}{g}^4\int d^{4}x{e}^{-iqx}{\int }^{\tilde{A}(t_f)=A(t_f)}D{\tilde{A}}_{\mu }D{A}_{\mu }{\int }^{\tilde{\psi }(t_f)=\psi (t_f)}D\widetilde{\stackrel{ˉ}{\psi }}D\tilde{\psi }D\bar{\psi }D\psi {\Psi }_{p_A}^{\ast }(\overrightarrow{\stackrel{~}{A}}(t_i),\tilde{\psi }(t_i))$$

$$\times {\Psi }_{p_B}^{\ast }(\overrightarrow{\stackrel{~}{A}}(t_i),\tilde{\psi }(t_i))e^{-i{S}_{\mathrm{QCD}}(\tilde{A},\tilde{\psi })}{e}^{i{S}_{\mathrm{QCD}}(A,\psi )}{\tilde{F}}^2(x)F^2(0){\Psi }_{p_A}(\vec{A}(t_i),\psi (t_i)){\Psi }_{p_B}(\vec{A}(t_i),\psi (t_i))$$

$$W(p_A,p_B,q)=g^4\int d^{4}x{e}^{-iqx}{\int }^{\tilde{A}(t_f)=A(t_f)}D{\tilde{A}}_{\mu }D{A}_{\mu }{\int }^{{\tilde{\psi }}_a(t_f)={\psi }_a(t_f)}D{\bar{\psi }}_{a}D{\psi }_{a}D{\widetilde{\stackrel{ˉ}{\psi }}}_{a}D{\tilde{\psi }}_a$$

$$\times e^{-i{S}_{\mathrm{QCD}}(\tilde{A},{\tilde{\psi }}_a)}{e}^{i{S}_{\mathrm{QCD}}(A,{\psi }_a)}{\Psi }_{p_A}^{\ast }(\overrightarrow{\stackrel{~}{A}}(t_i),{\tilde{\psi }}_a(t_i)){\Psi }_{p_A}(\vec{A}(t_i),\psi (t_i))$$

$$\times {\int }^{\tilde{B}(t_f)=B(t_f)}D{\tilde{B}}_{\mu }D{B}_{\mu }{\int }^{{\tilde{\psi }}_b(t_f)={\psi }_b(t_f)}D{\bar{\psi }}_{b}D{\psi }_{b}D{\widetilde{\stackrel{ˉ}{\psi }}}_{b}D{\tilde{\psi }}_b$$

$$\times e^{-i{S}_{\mathrm{QCD}}(\tilde{B},{\tilde{\psi }}_b)}{e}^{i{S}_{\mathrm{QCD}}(B,{\psi }_b)}{\Psi }_{p_B}^{\ast }(\overrightarrow{\stackrel{~}{B}}(t_i),{\tilde{\psi }}_b(t_i)){\Psi }_{p_B}(\vec{B}(t_i),{\psi }_b(t_i))$$

$$\times \int D{C}_{\mu }{\int }^{\tilde{C}(t_f)=C(t_f)}D{\tilde{C}}_{\mu }\int D{\bar{\psi }}_{C}D{\psi }_C{\int }^{{\tilde{\psi }}_c(t_f)={\psi }_c(t_f)}D{\widetilde{\stackrel{ˉ}{\psi }}}_{C}D{\tilde{\psi }}_C{\tilde{F}}_C^2(x)F_C^2(0)e^{-i{\tilde{S}}_C+i{S}_C}$$

$$\int D{C}_{\mu }{\int }^{\tilde{C}(t_f)=C(t_f)}D{\tilde{C}}_{\mu }\int D{\bar{\psi }}_{C}D{\psi }_C{\int }^{{\tilde{\psi }}_c(t_f)={\psi }_c(t_f)}D{\widetilde{\stackrel{ˉ}{\psi }}}_{C}D{\tilde{\psi }}_C{g}^4{\tilde{F}}_C^2(x)F_C^2(0)e^{-i{\tilde{S}}_C+i{S}_C}$$

$$=e^{S_{\mathrm{eff}}(U,V,\tilde{U},\tilde{V})}\mathcal{O}(q,x;A,\tilde{A},{\psi }_a{\tilde{\psi }}_a;B,\tilde{B},{\psi }_b,{\tilde{\psi }}_b)$$

$$\displaystyle S_{\rm eff}(U,V,{\tilde{U}},{\tilde{V}})~{}=~{}2{\rm Tr}\int\!d^{2}x_{\perp}\big{[}-i{\tilde{U}}_{i}{\tilde{V}}^{i}+iU_{i}V^{i}$$

$$\displaystyle+~{}\big{(}\tilde{\cal L}_{i}({\tilde{U}},{\tilde{V}})\tilde{\cal L}^{i}({\tilde{U}},{\tilde{V}})-2\tilde{\cal L}_{i}({\tilde{U}},{\tilde{V}}){\cal L}^{i}(U,V)+{\cal L}^{i}(U,V){\cal L}^{i}(U,V)\big{)}\ln\sigma_{a}\sigma_{b}s+O(\ln\sigma_{a}\sigma_{b}s)^{2}\big{]}$$

$$U(x_{\perp })=[\infty p_2+x_{\perp },-\infty p_2+x_{\perp }{]}^{A_{\ast }},U_i=U^{†}i{\partial }_{i}U$$

$$V(x_{\perp })=[\infty p_1+x_{\perp },-\infty p_1+x_{\perp }{]}^{B_{\bullet }},V_i=V^{†}i{\partial }_{i}V$$

$$W(p_A,p_B,q)=\int d^{4}x{e}^{-iqx}{\int }^{\tilde{A}(t_f)=A(t_f)}D{\tilde{A}}_{\mu }D{A}_{\mu }{\int }^{{\tilde{\psi }}_a(t_f)={\psi }_a(t_f)}D{\bar{\psi }}_{a}D{\psi }_{a}D{\widetilde{\stackrel{ˉ}{\psi }}}_{a}D{\tilde{\psi }}_a$$

$$\times e^{-i{S}_{\mathrm{QCD}}(\tilde{A},{\tilde{\psi }}_a)}{e}^{i{S}_{\mathrm{QCD}}(A,{\psi }_a)}{\Psi }_{p_A}^{\ast }(\overrightarrow{\stackrel{~}{A}}(t_i),{\tilde{\psi }}_a(t_i)){\Psi }_{p_A}(\vec{A}(t_i),\psi (t_i))$$

$$\times {\int }^{\tilde{B}(t_f)=B(t_f)}D{\tilde{B}}_{\mu }D{B}_{\mu }{\int }^{{\tilde{\psi }}_b(t_f)={\psi }_b(t_f)}D{\bar{\psi }}_{b}D{\psi }_{b}D{\widetilde{\stackrel{ˉ}{\psi }}}_{b}D{\tilde{\psi }}_b$$

$$\times e^{-i{S}_{\mathrm{QCD}}(\tilde{B},{\tilde{\psi }}_b)}{e}^{i{S}_{\mathrm{QCD}}(B,{\psi }_b)}{\Psi }_{p_B}^{\ast }(\overrightarrow{\stackrel{~}{B}}(t_i),{\tilde{\psi }}_b(t_i)){\Psi }_{p_B}(\vec{B}(t_i),{\psi }_b(t_i))$$

$$\times e^{S_{\mathrm{eff}}(U,V,\tilde{U},\tilde{V})}\mathcal{O}(q,x;A,{\psi }_a,\tilde{A},{\tilde{\psi }}_a;B,{\psi }_b,\tilde{B},{\tilde{\psi }}_b)$$

$$\tilde{A}(t_f)=A(t_f),{\tilde{\psi }}_a(t_f)={\psi }_a(t_f)\mathrm{and}\tilde{B}(t_f)=B(t_f),{\tilde{\psi }}_b(t_f)={\psi }_b(t_f)$$

$$A(x_{\bullet },x_{\perp })=\tilde{A}(x_{\bullet },x_{\perp }),{\psi }_a(x_{\bullet },x_{\perp })={\tilde{\psi }}_a(x_{\bullet },x_{\perp })$$

$$\mathrm{and}$$

$$B(x_{\ast },x_{\perp })=\tilde{B}(x_{\ast },x_{\perp }),{\psi }_b(x_{\ast },x_{\perp })={\tilde{\psi }}_b(x_{\ast },x_{\perp })$$

$$\sum _X⟨p_A,p_B|X⟩⟨X|p_A,p_B⟩=1$$

$$⟨p_A,p_B|e^{S_{\mathrm{eff}}(U,V,\tilde{U},\tilde{V})}|p_A,p_B⟩=1$$

$$\int D{C}_{\mu }{\int }^{\tilde{C}(t_f)=C(t_f)}D{\tilde{C}}_{\mu }\int D{\bar{\psi }}_{C}D{\psi }_C{\int }^{{\tilde{\psi }}_c(t_f)={\psi }_c(t_f)}D{\widetilde{\stackrel{ˉ}{\psi }}}_{C}D{\tilde{\psi }}_C{g}^4{\tilde{F}}_C^2(x)F_C^2(0)e^{-i{\tilde{S}}_C+i{S}_C}$$

$$=\mathcal{O}(q,x;A,{\psi }_a;B,{\psi }_b)$$

$$W(p_A,p_B,q)=\int d^{4}x{e}^{-iqx}⟨p_A|⟨p_B|\hat{\mathcal{O}}(q,x;\hat{A},{\hat{\psi }}_a;\hat{B},{\hat{\psi }}_b)|p_A⟩|p_B⟩$$

$$\hat{\mathcal{O}}(q,x;\hat{A},{\hat{\psi }}_a;\hat{B},{\hat{\psi }}_b)=\sum _{m,n}\int d{z}_{m}d{z}_n^{\prime }{c}_{m,n}^{\mu \nu }(q,x){\hat{\Phi }}_A(z_m){\hat{\Phi }}_B(z_n^{\prime })$$

$$W=\int d^{4}x{e}^{-iqx}\sum _{m,n}\int d{z}_m{c}_{m,n}^{\mu \nu }(q,x)⟨p_A|{\hat{\Phi }}_A(z_m)|p_A⟩\int d{z}_n^{\prime }⟨p_B|{\hat{\Phi }}_B(z_n^{\prime })|p_B⟩$$

$$\frac{64/s^2}{N_c^2-1}\int d^{4}x{e}^{-iqx}⟨p_A|{\hat{U}}_{\ast }^{mi}(x_{\bullet },x_{\perp }){\hat{U}}_{\ast }^{mj}(0)|p_A⟩⟨p_B|{\hat{V}}_{\bullet i}^n(x_{\ast },x_{\perp }){\hat{V}}_{\bullet j}^n(0)|p_B⟩$$

$${\hat{U}}_{\ast i}^a(z_{\bullet },z_{\perp })\equiv [-{\infty }_{\bullet },z_{\bullet }{]}_z^{ab}g{\hat{F}}_{\ast i}^b(z_{\bullet },z_{\perp }),{\hat{V}}_{\bullet i}^a(z_{\ast },z_{\perp })\equiv [-{\infty }_{\ast },z_{\ast }{]}_z^{ab}g{\hat{F}}_{\bullet i}^b(z_{\ast },z_{\perp })$$

Full text

††preprint: BNL-113982-2017-JA, JLAB-THY-17-2484††institutetext: a Physics Department, Old Dominion University, Norfolk, VA 23529, USA and Thomas Jefferson National Accelerator Facility, Newport News, VA 23606, USA††institutetext: b Physics Department, Brookhaven National Laboratory, Upton, NY 11973, USA

1 Introduction

Particle production in hadron-hadron scattering with transverse momentum of produced particle much smaller than its invariant mass is described in the framework of TMD factorization Collins:2011zzd ; Collins:1981uw ; Collins:1984kg ; Ji:2004wu ; GarciaEchevarria:2011rb . The typical example is the Higgs production at LHC through gluon-gluon fusion. Factorization formula for particle production in hadron-hadron scattering looks like Collins:2011zzd ; Collins:2014jpa

[TABLE]

where $\eta$ is the rapidity, ${\cal D}_{f/A}(x,z_{\perp},\eta)$ is the TMD density of a parton $f$ in hadron $A$, and $\sigma(ff\rightarrow H)$ is the cross section of production of particle $H$ of invariant mass $m_{H}^{2}=Q^{2}$ in the scattering of two partons. (For simplicity, we consider the scattering of unpolarized hadrons.)

In this paper we calculate the first power corrections $\sim{q_{\perp}^{2}\over Q^{2}}$ in a sense that we represent them as a TMD-like matrix elements of higher-twist operators. It should be noted that our method works for arbitrary relation between $s$ and $Q^{2}$ and between $q_{\perp}^{2}$ and hadron mass $m^{2}$ (provided that pQCD is applicable), but in this paper we only present the result for the physically interesting region $s\gg Q^{2}\gg q_{\perp}^{2}\gg m^{2}$.

To obtain formula (1) with first corrections we use factorization in rapidity Balitsky:1998ya . We denote quarks and gluons with rapidity close to the rapidity of the projectile and target protons as $A$-fields and $B$-fields, respectively. We call the remaining fields in the central region of rapidity by the name $C$-fields and integrate over them in the corresponding functional integral. At this step, we get the effective action depending on $A$ and $B$ fields. The subsequent integration over $A$ fields gives matrix elements of some TMD-like operators switched between projectile proton states and integration over $B$ fields will give matrix elements between target states. 111It should be noted that due to the kinematics $Q^{2}\gg Q_{\perp}^{2},m^{2}$ we will not need the explicit form of the high-energy effective action which is much sought after in the small-x physics but not known up to now except a couple of first perturbative terms Balitsky:1998ya ; Balitsky:2005we ; Hatta:2005rn ; Hatta:2005ia ; Altinoluk:2009je .

The paper is organized as follows. In Sect. 2 we derive the TMD factorization from the double functional integral for the cross section of particle production. In Sect 2, which is central to our approach, we explain the method of calculation of higher-twist power corrections based on a solution of classical Yang-Mills equations. In Sect. 4 we find the leading higher-twist correction to particle production in the region $s\gg Q^{2}\gg q_{\perp}^{2}$. Finally, in Sect. 5 we compare our calculations in the small-$x$ limit to the classical field resulting from the scattering of two shock waves. The Appendices contain proofs of some necessary technical statements.

2 TMD factorization from functional integral

We consider production of an (imaginary) scalar particle $\Phi$ in proton-proton scattering. This particle is connected to gluons by the vertex

[TABLE]

This is a ${m_{H}\over m_{t}}\ll 1$ approximation Kniehl:1995tn ; Shifman:1979eb for Higgs production via gluon fusion at LHC with

[TABLE]

where $\alpha_{s}={g^{2}\over 4\pi}$ as usual. 222For finite $m_{t}$ the constant $g_{H}$ should be multiplied by {3\tau\over 2}\big{[}1+(1-\tau)\arcsin^{2}{1\over\sqrt{\tau}}\big{]} with $\tau={4m_{t}^{2}\over m_{H}^{2}}$ Gunion:1989we . The differential cross section of $\Phi$ production has the form

[TABLE]

where we defined the “hadronic tensor” $W(p_{A},p_{B},q)$ as

[TABLE]

As usual, $\sum_{X}$ denotes the sum over full set of “out” states. It can be represented by double functional integral

[TABLE]

Here the fields $A,\psi$ correspond to the amplitude $\langle X|F^{2}(0)|p_{A},p_{B}\rangle$, fields ${\tilde{A}},\tilde{\psi}$ correspond to complex conjugate amplitude $\langle p_{A},p_{B}|F^{2}(x)|X\rangle$ and $\Psi_{p}(\vec{A}(t_{i}),\psi(t_{i}))$ denote the proton wave function at the initial time $t_{i}$. The boundary conditions ${\tilde{A}}(t_{f})=A(t_{f})$ and $\tilde{\psi}(t_{f})=\psi(t_{f})$ reflect the sum over all states $X$, cf. Refs. Balitsky:1988fi , Balitsky:1990ck , Balitsky:1991yz .

We use Sudakov variables $p=\alpha p_{1}+\beta p_{2}+p_{\perp}$ and the notations $x_{\bullet}\equiv x_{\mu}p_{1}^{\mu}$ and $x_{\ast}\equiv x_{\mu}p_{2}^{\mu}$ for the dimensionless light-cone coordinates ($x_{\ast}=\sqrt{s\over 2}x_{+}$ and $x_{\bullet}=\sqrt{s\over 2}x_{-}$). Our metric is $g^{\mu\nu}~{}=~{}(1,-1,-1,-1)$ so that $p\cdot q~{}=~{}(\alpha_{p}\beta_{q}+\alpha_{q}\beta_{p}){s\over 2}-(p,q)_{\perp}$ where $(p,q)_{\perp}\equiv-p_{i}q^{i}$. Throughout the paper, the sum over the Latin indices $i$, $j$… runs over the two transverse components while the sum over Greek indices runs over the four components as usual.

To derive the factorization formula, we separate the (quark and gluon) fields in the functional integral (5) into three sectors: “projectile” fields $A_{\mu},\psi_{a}$ with $|\beta|<\sigma_{a}$, “ target” fields with $|\alpha|<\sigma_{b}$ and “central rapidity” fields $C_{\mu},\psi$ with $|\alpha|>\sigma_{b}$ and $|\beta|>\sigma_{a}$: 333The standard factorization scheme for particle production in hadron-hadron scattering is splitting the diagrams in collinear to projectile part, collinear to target part, hard factor, and soft factor Collins:2011zzd . Here we factorize only in rapidity. For our purpose of calculation of power corrections in the tree approximation it is sufficient; however, we hope to treat possible logs of transverse scales in loop corrections in the same way as it was done in our rapidity evolution equations for gluon TMDs in Refs. Balitsky:2015qba ; Balitsky:2016dgz .

[TABLE]

where $S_{C}=S_{\rm QCD}(A+B+C)-S_{\rm QCD}(A)-S_{\rm QCD}(B)$.

Our goal is to integrate over central fields and get the amplitude in the factorized form, as a (sum of) products of functional integrals over $A$ fields representing projectile matrix elements (TMDs) and functional integrals over $B$ fields representing target matrix elements. In the spirit of background-field method, we “freeze” projectile and target fields (and denote them the ${\bar{A}}$, ${\bar{\xi}}_{a}$, $\xi_{a}$ and ${\bar{B}}$, ${\bar{\xi}}_{b}$, $\xi_{b}$ respectively) and get a sum of diagrams in these external fields. Since $|\beta|<\sigma_{a}$ in the projectile fields and $|\alpha|<\sigma_{b}$ in the target fields, at the tree-level one can set with power accuracy $\beta=0$ for the projectile fields and $\alpha=0$ for the target fields - the corrections will be O\big{(}{m^{2}\over\sigma_{a}s}\big{)} and O\big{(}{m^{2}\over\sigma_{b}s}\big{)}. Beyond the tree level, one should expect that the integration over $C$ fields will produce the logarithms of the cutoffs $\sigma_{a}$ and $\sigma_{b}$ which will cancel with the corresponding logs in gluon TMDs of the projectile and the target.

As usual, diagrams disconnected from the vertices $F^{2}(x)$ and $F^{2}(0)$ (“vacuum bubbles” in external fields) exponentiate so the result has the schematic form

[TABLE]

where ${\cal O}^{\mu\nu}(q,x;A,\psi_{A};B,\psi_{B})$ is a sum of diagrams connected to ${\tilde{F}}^{2}(x)F^{2}(0)$. Since rapidities of central fields and $A$, $B$ fields are very different, one should expect the result of integration over C-fields to be represented in terms of Wilson-line operators constructed form $A$ and $B$ fields.

The effective action has the form

[TABLE]

where Wilson lines $U$ are made from projectile fields

[TABLE]

and Wilson lines $V$ from target fields

[TABLE]

and similarly for ${\tilde{U}}$ and ${\tilde{V}}$ in the left sector. The explicit form of “Lipatov vertices” $L_{i}(U,V)$ is presented in Balitsky:2004rr . Unfortunately, the effective action beyond the first two terms in (8) is unknown, but we will demonstrate below that for our purposes we do not need the explicit form of the effective action.

After integration over $C$ fields the amplitude (5) can be rewritten as

[TABLE]

Note that due to boundary conditions at $t_{f}$ in the above integral, the functional integral over $C$ fields in Eq. (7) should be done in the background of the $A$ and $B$ fields satisfying

[TABLE]

Our approximation at the tree level is that $\beta=0$ for $A,{\tilde{A}}$ fields and $\alpha=0$ for $B,{\tilde{B}}$ fields which corresponds to $A=A(x_{\bullet},x_{\perp}),~{}{\tilde{A}}={\tilde{A}}(x_{\bullet},x_{\perp})$ and $B=B(x_{\ast},x_{\perp}),~{}{\tilde{B}}={\tilde{B}}(x_{\ast},x_{\perp})$.

Now comes the important point: because of boundary conditions (10), for the purpose of calculating the integral (7) over central fields one can set

[TABLE]

Indeed, because $A,\psi$ and ${\tilde{A}},\tilde{\psi}$ do not depend on $x_{\ast}$, if they coincide at $x_{\ast}=\infty$ they should coincide everywhere. Similarly, if $B,\psi_{b}$ and ${\tilde{B}},\tilde{\psi}_{b}$ do not depend on $x_{\bullet}$, if they coincide at $x_{\bullet}=\infty$ they should be equal.

It should be emphasized that the boundary conditions (10) mean the summation over all intermediate states in corresponding projectile and target matrix elements in the functional integrals over projectile and target fields. Without the sum over all intermediate states the conditions (11) are no longer true. For example, if we would like to measure another particle or jet in the fragmentation region of the projectile, the second condition in Eq. (11) breaks down.

Next important observation is that due to Eqs. (11) the effective action (8) vanishes for background fields satisfying conditions (10). For the first two terms displayed in (8) it is evident, but it is easy to see that the effective action in the background fields satisfying (11) should vanish due to unitarity. Indeed, let us consider the functional integral (5) without sources $\tilde{F}^{2}(x)F^{2}(0)$. It describes the matrix element (12) without $\Phi$ production, that is

[TABLE]

(modulo appropriate normalization of $|p_{A}\rangle$ and $|p_{B}\rangle$ states). If we perform the same decomposition into $A$, $B$, and $C$ fields as in Eq. (5) we will see integral (9) without ${\cal O}^{\mu\nu}(q,x,y;A,\psi_{a},{\tilde{A}},\tilde{\psi}_{a};B,\psi_{b},{\tilde{B}},\tilde{\psi}_{b})$ which can be represented as

[TABLE]

which means that the effective action should vanish for the Wilson-line operators constructed from the fields satisfying Eqs. (11). Summarizing, we see that at the tree level in our approximation

[TABLE]

where now $S_{C}~{}=~{}S_{\rm QCD}(C+A+B)-S_{\rm QCD}(A)-S_{\rm QCD}(B)$ and ${\tilde{S}}_{C}~{}=~{}S_{\rm QCD}({\tilde{C}}+A+B)-S_{\rm QCD}(A)-S_{\rm QCD}(B)$. It is known that in the tree approximation the double functional integral (14) is given by a set of retarded Green functions in the background fields Gelis:2003vh ; Gelis:2006yv ; Gelis:2007kn (see also Appendix A for the proof). Since the double functional integral (14) is given by a set of retarded Green functions (in the background field $A+B$), the calculation of tree-level contributions to, say, $F^{2}(x)$ in the r.h.s. of Eq. (14) is equivalent to solving YM equation for $A_{\mu}(x)$ (and $\psi(x)$) with boundary conditions that the solution has the same asymptotics at $t\rightarrow-\infty$ as the superposition of incoming projectile and target background fields.

The hadronic tensor (9) can now be represented as

[TABLE]

where $\hat{\cal O}(q,x;{\hat{A}},{\hat{\psi}}_{a};{\hat{B}},{\hat{\psi}}_{b})$ should be expanded in a series in ${\hat{A}},{\hat{\psi}}_{a};{\hat{B}},{\hat{\psi}}_{b}$ operators and evaluated between the corresponding (projectile or target) states: if

[TABLE]

(where $c^{\mu\nu}_{m,n}$ are coefficients and $\Phi$ can be any of $A_{\mu}$, $\psi$ or ${\bar{\psi}}$) then 444 Our logic here is the following: to get the expression for $\hat{\cal O}$ in Eq. (14) we calculate ${\cal O}$ in the background of two external fields $\Phi_{A}=(A_{\mu},\psi_{a})$ and $\Phi_{B}=(B_{\mu},\psi_{b})$ and then promote them to operators $\hat{\Phi}_{A}$ and $\hat{\Phi}_{B}$ in the obtained expressions for ${\cal O}$. However, there is a subtle point in the promotion of background fields to operators. When we are calculating ${\cal O}$ as the r.h.s. of Eq. (14) the fields $\Phi_{A}$ and $\Phi_{B}$ are c-numbers; on the other hand, after functional integration in Eq. (5) they become operators which must be time-ordered in the right sector and anti-time-ordered in the left sector. Fortunately, as we shall see below, all these operators are separated either by space-like distances or light-cone distances so all of them (anti) commute and thus can be treated as $c$-numbers.

[TABLE]

As we will demonstrate below, the relevant operators are quark and gluon fields with Wilson-line type gauge links collinear to either $p_{2}$ for $A$ fields or $p_{1}$ for $B$ fields.

3 Power corrections and solution of classical YM equations

3.1 Power counting for background fields

As we discussed in previous Section, to get the hadronic tensor in the form (15) we need to calculate the functional integral (14) in the background of the fields (11). To understand the relative strength of Lorentz components of these fields, let us compare the typical term in the leading contribution to $W$

[TABLE]

where

[TABLE]

and some typical higher-twist terms. As we mentioned, we consider $W(p_{A},p_{B},q)$ in the region where $s,Q^{2}\gg Q_{\perp}^{2},m^{2}$ while the relation between $Q_{\perp}^{2}$ and $m^{2}$ and between $Q^{2}$ and $s$ may be arbitrary. So, for the purpose of counting of powers of $s$, we will not distinguish between $s$ and $Q^{2}$ (although at the final step we will be able to tell the difference since our final expressions for higher-twist corrections will have either $s$ or $Q^{2}$ in denominators). Similarly, for the purpose of power counting we will not distinguish between $m$ and $Q_{\perp}$ and will introduce $m_{\perp}$ which may be of order of $m$ or $Q_{\perp}$ depending on matrix element.

The estimate of the leading-twist matrix element between projectile states is

[TABLE]

(here we assume normalization $\langle p_{A}|p_{A}\rangle~{}=~{}1$ for simplicity).

The typical higher-twist correction is proportional to (see e.g. Eq. (64))

[TABLE]

Since $x^{\perp}_{i}\sim{q^{\perp}_{i}\over q_{\perp}^{2}}\sim{1\over m_{\perp}}$ we see that an extra ${\hat{F}}_{\mu i}$ in the matrix element between projectile states brings $p_{1\mu}m_{\perp}$ which means that ${\hat{U}}_{\ast i}\sim sm_{\perp}$.

Next, some of the higher-twist matrix elements have an extra $U_{kl}$ like

[TABLE]

where

[TABLE]

Since we consider only unpolarized projectile and target hadrons

[TABLE]

and, comparing this to Eq. (20), we see that an extra ${\hat{F}}_{kl}$ can bring an extra $m_{\perp}^{2}$. Combining this with an estimate $U_{\ast i}\sim sm_{\perp}$ we see that the typical field ${\bar{A}}_{\ast}$ is of order $s$ while ${\bar{A}}_{i}\sim m_{\perp}$. Similarly, for the target fields we get ${\bar{B}}_{\bullet}\sim s$, ${\bar{B}}_{i}\sim m_{\perp}$.

Some of the power corrections involve matrix elements like

[TABLE]

where

[TABLE]

An extra field strength operator ${\hat{F}}^{\mu\nu}$ between the projectile states can bring ${p_{A}^{\mu}p_{2}^{\nu}\over p_{A}\cdot p_{2}}-\mu\leftrightarrow\nu$ so that ${\hat{F}}_{\ast{\bullet}}\sim sm^{2}$ 555The denominator $p_{A}\cdot p_{2}$ is due to the fact that $p_{2}$ enters only through the direction of Wilson line and therefore the matrix element should not change under rescaling $p_{2}\rightarrow\lambda p_{2}$. Since ${\bar{A}}_{\ast}\sim s$ we see that ${\bar{A}}_{\bullet}\sim m^{2}_{\perp}$. Similarly, for the target we get ${\bar{B}}_{\ast}\sim m_{\perp}^{2}$.

Summarizing, the relative strength of the background gluon fields in projectile and target is

[TABLE]

To finish power counting, we need also the relative strength of quark background fields $\psi_{a}$ and $\psi_{b}$. From classical equations for projectile and target

[TABLE]

we get

[TABLE]

Thus, to find TMD factorization at the tree level (with higher-twist corrections) we need to calculate the functional integral (5) in the background fields of the strength given by Eqs. (27) and (29).

3.2 Approximate solution of classical equations

As we discussed in Sect 2, the calculation of the functional integral (14) over $C$-fields in the tree approximation reduces to finding fields $C_{\mu}$ and $\psi_{c}$ as solutions of Yang-Mills equations for the action $S_{C}~{}=~{}S_{\rm QCD}(C+A+B)-S_{\rm QCD}(A)-S_{\rm QCD}(B)$

[TABLE]

As we discussed above (see also Appendix A) the solution of Eq. (30) which we need corresponds to the sum of set of diagrams in background field ${\bar{A}}+{\bar{B}}$ with retarded Green functions (see Fig. 3).

The retarded Green functions (in the background-Feynman gauge) are defined as

[TABLE]

where

[TABLE]

and similarly for quarks.

The solutions of Eqs. (30) in terms of retarded Green functions give fields $C_{\mu}$ and $\psi_{c}$ that vanish at $t\rightarrow-\infty$. Thus, we are solving the usual classical YM equations

[TABLE]

with boundary conditions

[TABLE]

following from $C_{\mu},\psi_{c}\stackrel{{\scriptstyle t\rightarrow-\infty}}{{\rightarrow}}0$. These boundary conditions reflect the fact that at $t\rightarrow-\infty$ we have only incoming hadrons with “A” and “B” fields.

The solution of YM equations (33) in general case is yet unsolved problem, especially important for scattering of two heavy nuclei in semiclassical approximation. Fortunately, for our case of particle production with ${q_{\perp}\over Q}\ll 1$ we can construct the approximate solution of (33) as a series in this small parameter. However, before doing this, it is convenient to perform a gauge transformation so that the incoming projectile and target fields will no longer have large components $\sim s$ as ${\bar{A}}_{\ast}$ and ${\bar{B}}_{\bullet}$ in Eq. (27). Let us perform the gauge transformation of Eq. (33) and initial conditions (34) with the gauge matrix $\Omega(x)$ such that

[TABLE]

The existence of such matrix is proved in Appendix B by explicit construction. After such gauge transformation, the YM equation of course stays the same but the initial conditions (34) turn to

[TABLE]

where

[TABLE]

and $\Sigma_{a},\Sigma_{b}$ are defined as

[TABLE]

The initial conditions (36) look like the projectile fields in the light-like gauge $p_{2}^{\mu}A_{\mu}=0$ and target fields in the light-like gauge $p_{1}^{\mu}A_{\mu}=0$ so our construction of matrix $\Omega$ in a way proves that we can take the sum of projectile fields in one gauge and target fields in another gauge as a zero-order approximation for iterative solution of the YM equations. Note also that our power counting discussed in previous Section means that

[TABLE]

so we do not have large background fields $\sim s$ after this gauge transformation. Finally, the classical equations for projectile and target fields in this gauge read 666Here we consider only $u$, $d$, and $s$ quarks which can be regarded as massless.:

[TABLE]

where $U_{\mu\nu}\equiv\partial_{\mu}U_{\nu}-\partial_{\nu}U_{\mu}-i[U_{\mu},U_{\nu}]$, $D_{U}^{\mu}\equiv(\partial^{\mu}-i[U^{\mu},)$ and similarly for $V$ fields.

We will solve Eqs. (33) iteratively, order by order in perturbation theory, starting from the zero-order approximation in the form of the sum of projectile and target fields

[TABLE]

and improving it by calculation of Feynman diagrams with retarded propagators in the background fields (41).

The first step is the calculation of the linear term for the trial configuration (41). We rewrite field strength components as

[TABLE]

Note that $U_{\ast i}\sim V_{{\bullet}i}\sim sm_{\perp}$, $U_{\ast{\bullet}}\sim V_{\ast{\bullet}}\sim sm_{\perp}^{2}$ while all other components are not large.

The linear term has the form

[TABLE]

where ${\cal D}^{j}~{}\equiv~{}\partial^{j}-iU^{j}-iV^{j}$, ${\cal D}_{\bullet}=\partial_{\bullet}-iU_{\bullet}$, and ${\cal D}_{\ast}=\partial_{\ast}-iV_{\ast}$. The power-counting estimates for linear terms in Eq. (43) are

[TABLE]

Note that the order of perturbation theory is labeled by $(...)^{[n]}$ and the order of expansion in the parameter ${m_{\perp}^{2}\over s}$ by $(...)^{(n)}$.

With the linear term (43), a couple of first terms in perturbative series are

[TABLE]

for gluon fields (in the background-Feynman gauge) and

[TABLE]

for quarks where

[TABLE]

are operators in external zero-order fields (41). Hereafter we use Schwinger’s notations for propagators in external fields normalized according to (x|F(p)|y)\equiv\int\!{\textstyle d}\lower 0.12915pt\hbox{\kern-3.80005pt{}^{\scriptstyle-}}\kern-0.50003pt{}^{4}p~{}e^{-ip(x-y)}F(p) . Moreover, when it will not lead to a confusion, we will use short-hand notation

${1\over{\cal O}}{\cal O}^{\prime}(x)~{}\equiv~{}\int\!d^{4}z~{}(x|{1\over{\cal O}}|z){\cal O}^{\prime}(z)$. Next iterations will give us a set of tree-level Feynman diagrams in the background field $U_{\mu}+V_{\mu}$ and $\Sigma_{a}+\Sigma_{b}$.

Let us consider the fields in the first order in perturbation theory:

[TABLE]

Here $\alpha$, $\beta$, and $p_{\perp}$ are understood as differential operators $\alpha=i{\partial\over\partial x_{\bullet}}$, $\beta=i{\partial\over\partial x_{\ast}}$ and $p_{i}=i{\partial\over\partial x^{i}}$.

Now comes the central point of our approach. Let us expand quark and gluon propagators in powers of background fields, then we get a set of diagrams shown in Fig. 3. The typical bare gluon propagator in Fig. 3 is

[TABLE]

Since we do not consider loops of $C$-fields in this paper, the transverse momenta in tree diagrams are determined by further integration over projectile (“A”) and target (“B”) fields in Eq. (9) which converge on either $q_{\perp}$ or $m$. On the other hand, the integrals over $\alpha$ converge on either $\alpha_{q}$ or $\alpha\sim 1$ and similarly the characteristic $\beta$’s are either $\beta_{q}$ or $\sim 1$. Since $\alpha_{q}\beta_{q}s=Q_{\parallel}^{2}\gg Q_{\perp}^{2}$, one can expand gluon and quark propagators in powers of ${p_{\perp}^{2}\over\alpha\beta s}$

[TABLE]

The explicit form of operators ${1\over\alpha+i\epsilon}$, ${1\over\beta+i\epsilon}$, and ${1\over(\alpha+i\epsilon)(\beta+i\epsilon)}$ is

[TABLE]

After the expansion (50), the dynamics in the transverse space effectively becomes trivial: all background fields stand either at $x$ or at [math]. (This validates the reasoning in the footnote on page 3).

One may wonder why we do not cut the integrals in Eq. (51) to $|\alpha|>\sigma_{b}$ and $|\beta|>\sigma_{a}$ according to the definition of $C$ fields in Sect. 2. 777 Such cutoffs for integrals over $C$ fields are introduced explicitly in the framework of soft-collinear effecive theory, see the review Rothstein:2016bsq . The reason is that in the diagrams like Fig. 3 with retarded propagators (51) one can shift the contour of integration over $\alpha$ and/or $\beta$ to the complex plane away to avoid the region of small $\alpha$ or $\beta$. 888This may be wrong if there is pinching of poles in the integrals over $\alpha$ or $\beta$ but we will see that in our integrals for the tree-level power corrections the pinching of poles never occurs. In the higher orders in perturbation theory the pinching does occur so one needs to formulate a subtraction program to avoid double counting.

Note that the background fields are also smaller than typical $p_{\parallel}^{2}\sim s$. Indeed, from Eq. (39) we see that $p_{\bullet}={s\over 2}\beta\gg U_{\bullet}\sim m^{2}$ ( because $\alpha\geq\alpha_{q}\gg{m^{2}\over s}$) and similarly $p_{\ast}\gg V_{\ast}$. Also $(p_{i}+U_{i}+V_{i})^{2}~{}\sim~{}q_{\perp}^{2}\ll p_{\parallel}^{2}$. The only exception is the fields $V_{{\bullet}i}$ or $U_{\ast i}$ which are of order of $sm_{\perp}$ but we will see that effectively the expansion in powers of these fields is cut at the second term with our accuracy.

3.3 Twist expansion of classical gluon fields

Now we expand the classical gluon fields in powers of ${p_{\perp}^{2}\over p_{\parallel}^{2}}\sim{m_{\perp}^{2}\over s}$. It is clear that for the leading higher-twist correction we need to take into account only the first two terms (45) of the perturbative expansion of classical field. The expansion (45) of gluon field $A_{\bullet}$ takes the form

[TABLE]

where

[TABLE]

Similarly, from Eq. (45) one obtains

[TABLE]

and

[TABLE]

where ($n=1,2$)

[TABLE]

In these formulas the singularity in ${1\over\alpha}$ is always causal ${1\over\alpha+i\epsilon}$ and similarly for ${1\over\beta}\equiv{1\over\beta+i\epsilon}$ and ${1\over p_{\parallel}^{2}}\equiv{1/s\over(\alpha+i\epsilon)(\beta+i\epsilon)}$, see Eq. (51).

The corresponding expansion of field strengths reads

[TABLE]

Power corrections to hadronic tensor are proportional to

[TABLE]

so

[TABLE]

and the leading higher-twist correction is proportional to

[TABLE]

4 Leading higher-twist correction at $s\gg Q^{2}\gg Q_{\perp}^{2}\gg m^{2}$

As we mentioned in the Introduction, our method is relevant for calculation of higher-twist corrections at any $s,Q^{2}\gg Q_{\perp}^{2},m^{2}$. However, the expressions become manageable in the physically interesting case $s\gg Q^{2}\gg Q_{\perp}^{2}\gg m^{2}$ which we consider in this Section. 999We also assume that the scalar particle is emitted in the central region of rapidity so $\alpha_{q}s\sim\beta_{q}s\gg Q^{2}$. We will demonstrate that the leading correction in this region comes from the following part of Eq. (58)

[TABLE]

where

[TABLE]

The higher-twist correction coming from the second term in the r.h.s. will be $\sim{Q_{\perp}^{2}\over Q^{2}}$ whereas other terms in the r.h.s. of Eq. (58) yield contributions $\sim{Q_{\perp}^{2}\over s}$, $\sim{Q_{\perp}^{2}\over\alpha_{q}s}$, or $\sim{Q_{\perp}^{2}\over\beta_{q}s}$ all of which are small (see the footnote 9). In this approximation we get

[TABLE]

where the first term is the leading order and the second is the higher-twist correction.

Substituting our approximation (61) to Eq. (4) and promoting background fields to operators as discussed in Sect. 2 we get (note that $\alpha_{q}\beta_{q}s=Q_{\parallel}^{2}\simeq Q^{2}$):

[TABLE]

where we used formula Kaplan_Resnikoff_1967 ; MacFarlane:1968vc

[TABLE]

Since an extra $U_{\ast k}$ (or $V_{{\bullet}k}$) brings $s{x_{i}\over x_{\perp}^{2}}$ 101010To see this, we compared matrix elements of leading-twist operator $\langle p_{A}|U_{\ast}^{mi}(x_{\bullet},x_{\perp})U_{\ast}^{mj}(0)|p_{A}\rangle$ and higher-twist operator $\langle p_{A}|U^{a}_{\ast i}(x_{\bullet},x_{\perp})U^{b}_{\ast j}(x^{\prime}_{\bullet},x_{\perp})U_{\ast r}^{c}(0)|p_{A}\rangle$ between quark states which gives an extra $s{x_{r}\over x_{\perp}^{2}}$ modulo some logarithms.

we see that the higher-twist correction in the r.h.s of Eq. (64) is $\sim{q_{\perp}^{2}\over Q^{2}}$ so it gives the leading power correction in the region $s\gg Q^{2}=m_{\Phi}^{2}\gg q_{\perp}^{2}\gg m^{2}$. The TMD factorization formula with the higher-twist correction (64) is the main result of the present paper.

We parametrize gluon TMD for unpolarized protons as (cf. Ref. Mulders:2000sh )

[TABLE]

where $\sigma_{b}$ is the cutoff in $\alpha$ integration in the target matrix elements, see the discussion in Ref. Balitsky:2015qba . The normalization here is such that $D_{g}(\beta_{q},k_{\perp}^{2};\sigma_{b})$ is an unintegrated gluon distribution:

[TABLE]

where $D_{g}(\beta_{q},\mu^{2})$ is the usual gluon parton density (this formula is correct in the leading log approximation, see the discussion in Ref. Balitsky:2015qba ).

Next, the three-gluon matrix element in Eq. (64) for unpolarized hadrons can be parametrized as

[TABLE]

At large $k_{\perp}^{2}$ gluon TMDs in the r.h.s. of Eq. (66) behave as $D_{g}(\beta_{q},k^{2}_{\perp})\sim{1\over k_{\perp}^{2}}$ and $H_{g}(\beta_{q},k_{\perp}^{2})\sim{1\over k_{\perp}^{4}}$. Similarly, one should expect that $D_{i}^{g}(\beta_{q},k^{2}_{\perp})\sim{1\over k_{\perp}^{2}}$ and $H_{1}^{g}(\beta_{q},k_{\perp}^{2})\sim{1\over k_{\perp}^{4}}$.

It is well known that in our kinematic region $s\gg Q^{2}\gg Q_{\perp}^{2}$ gluon TMDs (66) possess Sudakov logs of the type

[TABLE]

One should expect double-logs of this type in $D_{i}^{g}(\beta_{q},k_{\perp}^{2};\sigma_{b})$ and $H_{1}^{g}(\beta_{q},k_{\perp}^{2};\sigma_{b})$, too.

Let us now demonstrate that the terms in $(F^{2}(x))^{(0)}$ (see Eq. (59)) which we neglected give small contributions. For example, consider the following contribution to $F^{2}(x)F^{2}(0)$:

[TABLE]

The corresponding contribution to hadronic tensor $W$ has the form

[TABLE]

Note that unlike Eq. (64), the factor in the denominator is $\alpha_{q}s\gg Q^{2}$ so the contribution (71) is power suppressed in comparison to Eq. (64) in our kinematic region. 111111Of course, this power suppression may be moderated by difference in logarithmic evolution of operators in the r.h.s.’s of Eqs. (64) and (71), but one should expect the evolution of these operators to be of the same order of magnitude.

As a less trivial example, consider the following term in $F^{2}(x)F^{2}(0)$

[TABLE]

The corresponding contribution to hadronic tensor $W$ reads

[TABLE]

where we used

[TABLE]

In both examples (71) and (73) the factor ${1\over\alpha_{q}}$ comes from an extra integration over $x^{\prime}_{\bullet}$ in $U_{i}$, see Eq. (37):

[TABLE]

The way to figure out such integrations is very simple: take $\alpha_{q}\rightarrow 0$ and check if there is an infinite integration of the type $\int_{-\infty}^{x_{\bullet}}dx^{\prime}_{\bullet}$. Evidently, it may happen if we have a single $U_{i}(x)$ (without any additional $U$-operators) at the point $x$, or a single $U_{i}(0)$.

Similarly, the factor ${1\over\beta_{q}}$ comes from an extra integration over $x^{\prime}_{\ast}$ in $V_{i}$ in Eq. (37) so an indication of such contribution is the infinite integration $\int_{-\infty}^{x_{\ast}}dx^{\prime}_{\ast}$ in the limit $\beta_{q}\rightarrow 0$ which translates to the condition of a single $V_{i}$ at the point $x$ or at the point [math].

Thus, to get the terms $\sim{1\over Q^{2}}$ we need to find contributions which satisfy both of the above conditions which singles out the contribution (63).

5 Small-x limit and scattering of shock waves

Let us consider the hadronic tensor

[TABLE]

in the small-x limit $s\rightarrow\infty$, $Q^{2}$ and $q_{\perp}^{2}$ - fixed. At first, let us not impose the condition $Q^{2}\gg q_{\perp}^{2}$ which means that the relation between $x_{\parallel}^{2}$ and $x_{\perp}^{2}$ is arbitrary (later we will see that $Q^{2}\gg q_{\perp}^{2}$ corresponds to $x_{\parallel}^{2}\ll x_{\perp}^{2}$).

The small-x limit may be obtained by rescaling $s\rightarrow\lambda^{2}s\Leftrightarrow p_{1}\rightarrow\lambda p_{1},p_{2}\rightarrow\lambda p_{2}$. As discussed in Refs. Balitsky:1995ub ; Balitsky:1998ya ; Balitsky:2004rr , the only components of field strength surviving in this rescaling are $U_{\ast i}(x_{\bullet},x_{\perp})$ and $V_{{\bullet}i}(x_{\ast},x_{\perp})$. Moreover, if we study classical fields at longitudinal distances which does not scale with $\lambda$, we can replace the projectile and target fields by infinitely thin “shock waves”

[TABLE]

However, since we need to compare the classical fields in the small-$x$ limit to our expressions (57) at small longitudinal distances, we will keep $x_{\ast}$ and $x_{\bullet}$ dependence for a while.

As described above, to find the classical fields we can start with the trial configuration

[TABLE]

with the linear term

[TABLE]

and improve it order by order in $L_{\mu}$. In this way we’ll get a set of Feynman diagrams in the background field (77). Unfortunately, in the general case of arbitrary relation between $q_{\parallel}$ and $q_{\perp}$ we no longer have a small parameter ${p_{\perp}^{2}\over p_{\parallel}^{2}}$ so we need explicit expressions for propagators in the background fields, and, in addition, we need all orders in the expansion of linear term (78). Still, we can compare our calculations with the perturbative expansion of classical fields in powers of the “parameter” $[U_{i},V_{j}]$ carried out in Refs. Balitsky:1998ya ; Balitsky:2005we . In the leading order in perturbation theory only the first line of Eq. (45) survives and we get

[TABLE]

The corresponding expressions for field strengths are

[TABLE]

In the last line we used the identity

[TABLE]

and the fact that in the small-$x$ limit $\partial_{i}U_{j}-\partial_{j}U_{i}-i[U_{i},U_{j}]~{}=~{}\partial_{i}V_{j}-\partial_{j}V_{i}-i[V_{i},V_{j}]~{}=~{}0$.

Let us discuss now how our approximation ${p_{\perp}^{2}\over p_{\parallel}^{2}}\ll 1$ looks in the coordinate space. The explicit expressions for fields (80) are

[TABLE]

where

[TABLE]

At longitudinal distances $x_{\bullet},x_{\ast}\sim 1$ these expressions agree with Eq. (52) from Ref. Balitsky:1998ya after the replacement (76).

Now let us compare the fields (82) at small longitudinal distances to our approximate solution (57). Let us start with $F_{ij}(x)$ in the last line in Eq. (82). If $(x-z)_{\parallel}^{2}$ is smaller than the characteristic transverse distances in the integral over $z_{\perp}$ one can replace $[U_{\ast i}(z_{\bullet},z_{\perp}),V_{{\bullet}j}(z_{\ast},z_{\perp})]$ by $[U_{\ast i}(z_{\bullet},x_{\perp}),V_{{\bullet}j}(z_{\ast},x_{\perp})]$ and get

[TABLE]

which is exactly the last line in Eq. (57). Similarly, the third line in Eq. (82) reproduces $F_{\ast{\bullet}}$ in the fourth line in Eq. (57).

Next, $gF^{(0)a}_{{\bullet}i}$ in second line in Eq. (57) in the leading order in perturbation theory turns to

[TABLE]

On the other hand, the first line in Eq. (82) at small $(x-z)_{\parallel}$ gives

[TABLE]

which agrees with Eq. (85) after integration over $z_{\perp}$. Similarly, one can check the consistency of two expressions for $F_{\ast i}$.

6 Conclusions and outlook

We have formulated the approach to TMD factorization based on the factorization in rapidity and found the leading higher-twist contribution to the production of a scalar particle (e.g. Higgs) by gluon-gluon fusion in the hadron-hadron scattering. Up to now our results are obtained in the tree-level approximation when the question of exact matching of cutoffs in rapidity does not arise. However, this question will become crucial starting from the first loop. In our previous papers we calculated the evolution of gluon TMD with respect to our rapidity cutoff so we need to match it to the coefficient functions in front of TMD operators. The work is in progress.

Also, we obtained power corrections for particle production only in the case of gluon-gluon fusion. It would be interesting (and we plan) to find power corrections to Drell-Yan process. There is a statement that for semi-inclusive deep inelastic scattering (SIDIS) the leading-order TMDs have different directions of Wilson lines: one to $+\infty$ and another to $-\infty$. We think that the same directions of Wilson lines will be in the case of power corrections and we plan to study this question in forthcoming publications.

The authors are grateful to J.C. Collins, S. Dawson, A. Kovner, D. Neill, A. Prokudin, T. Rogers, and R. Venugopalan for valuable discussions. This material is based upon work supported by the U.S. Department of Energy, Office of Science, Office of Nuclear Physics under contracts DE-AC02-98CH10886 and DE-AC05-06OR23177.

7 Appendix A

In this Section we prove that the field $C_{\mu}$ created by a source $J_{\mu}$ in the presence of external fields ${\bar{A}}_{\mu}$ and ${\bar{B}}_{\mu}$ 121212For simplicity, in this section we disregard quarks so in our case $J_{\mu}$ is Eq. (43) without quark terms.

[TABLE]

is given by a set of Feynman diagrams with retarded Green functions (note that Eq. (87) implies that $J_{\mu},{\bar{A}}_{\mu}$, and ${\bar{B}}_{\mu}$ are the same in the right and left part of the amplitude). Hereafter we use the notation $\Box_{\mu\nu}\equiv{\bar{P}}^{2}g_{\mu\nu}+2i{\bar{G}}_{\mu\nu}$.

First, we consider gluon propagators for the double functional integral over $C$ fields in the background filelds ${\bar{A}}=\bar{\tilde{A}}$, ${\bar{B}}=\bar{\tilde{B}}$ and prove that

[TABLE]

Note that we define $\langle\mathbb{O}\rangle$ in this Section as

[TABLE]

To prove Eq. (88), we write down

[TABLE]

and expand in powers of ${\cal O}_{\mu\nu}$.

In the trivial order Eqs. (88) immediately follow from the bare propagators for the double functional integral (89)

[TABLE]

where

[TABLE]

In the first order in ${\cal O}_{\mu\nu}$ we get

[TABLE]

so

[TABLE]

Similarly, it is easy to see that

[TABLE]

so

[TABLE]

In the second order in ${\cal O}_{\mu\nu}$

[TABLE]

so using the results (94) and (96) we get

[TABLE]

Similarly, it is easy to demonstrate that

[TABLE]

One can prove now Eq. (88) by induction using formulas

[TABLE]

Now we are in a position to prove Eq. (87). In the leading order in $g$ it is trivial: using Eqs. (88) one immediately sees that

[TABLE]

In the first order in $g$ (with one three-gluon vertex) we obtain

[TABLE]

which is the desired result.

Similarly, in the $g^{2}$ order one obtains after some algebra

[TABLE]

At arbitrary order in $g$ the structure similar to Eq. (103) can be proved by induction.

Thus, we see that Eq. (87) is given by a set of Feynman diagrams with retarded Green functions. In a similar way, one can demonstrate that

[TABLE]

8 Appendix B

To find matrix $\Omega(x)$ satisfying Eqs. (35) we will solve the following auxiliary problem: we fix $x_{\perp}$ as a parameter and find the solution of Yang-Mills equations

[TABLE]

in 2-dimensional gluodynamics with initial conditions

[TABLE]

Since 2-dimensional gluodynamics is a trivial theory, the solution of the equation (105) will be a pure-gauge field ${\cal A}_{\mu}=\Omega i\partial_{\mu}\Omega^{\dagger}$ with $\Omega(x_{\ast},x_{\bullet})$ being the sought-for matrix satisfying Eqs. (35).

Let us first demonstrate that the solution ${\cal A}_{\mu}(x_{\ast},x_{\bullet})$ of the YM equations (105) with boundary conditions (106) in two longitudinal dimensions is a pure gauge. To this end, we will construct ${\cal A}_{\mu}(x_{\ast},x_{\bullet})$ order by order in perturbation theory (see Fig. 3, but now in two dimensions) and prove that $F^{a}_{\mu\nu}({\cal A})~{}=~{}0$.

We are looking for the solution of Eq. (105) in the form

[TABLE]

Imposing the background-gauge condition

[TABLE]

we get the equation

[TABLE]

where ${\bar{D}}_{\mu}\equiv(\partial_{\mu}-ig[{\bar{A}}_{\mu},)$ and ${\bar{F}}_{\ast{\bullet}}~{}=~{}-ig[{\bar{A}}_{\ast},{\bar{B}}_{\bullet}]$. The boundary conditions (106) in terms of $C$ fields read

[TABLE]

It is convenient to rewrite the equation (109) in components as

[TABLE]

We will solve this equation by iterations in ${\bar{F}}_{\ast{\bullet}}$ and prove that ${\cal F}_{\ast{\bullet}}~{}=~{}0$ in all orders.

In the first order we get the equation

[TABLE]

The solution satisfying boundary conditions (110) has the form

[TABLE]

Using the explicit form of the propagators in external ${\bar{A}}_{\ast}$ and ${\bar{B}}_{\bullet}$ fields

[TABLE]

we get ${\bar{C}}^{(1)}$ in the form

[TABLE]

From this equation it is clear that $C_{\mu}^{(1)}(x_{\ast},x_{\bullet})$ vanishes if $x_{\ast}\rightarrow-\infty$ and/or $x_{\bullet}\rightarrow-\infty$ (recall that we assume ${\bar{A}}_{\ast}(x_{\bullet})\stackrel{{\scriptstyle x_{\bullet}\rightarrow\pm\infty}}{{\rightarrow}}0$ and ${\bar{B}}_{\bullet}(x_{\ast})\stackrel{{\scriptstyle x_{\ast}\rightarrow\pm\infty}}{{\rightarrow}}0$).

Also, form Eq. (113) we see that

[TABLE]

and therefore

[TABLE]

so in the first order in ${\bar{F}}$ the field strength of the solution of classical equation (109) vanishes.

In the second order the equations for the field $C_{\mu}$ take the form

[TABLE]

where we used Eq. (116) to reduce the r.h.s. Again, from the explicit form of the propagators (114) we get

[TABLE]

from which it is clear that ${\bar{C}}^{(2)}_{\mu}$ satisfy boundary conditions (110) (recall that we already proved that ${\bar{C}}^{(1)}_{\mu}$ satisfy Eq. (110)). Next, we use

[TABLE]

to prove that ${\cal F}_{\ast{\bullet}}$ vanishes in the second order:

[TABLE]

In the third order we get

[TABLE]

where again we used Eqs. (116) and (120) to reduce the r.h.s. The solution is

[TABLE]

Again, from the explicit form of propagators (114) it is clear that ${\bar{C}}_{\mu}^{(3)}$ satisfy boundary conditions (106) if ${\bar{C}}_{\mu}^{(1)}$ and ${\bar{C}}_{\mu}^{(2)}$ do (which we already proved). Next, from

[TABLE]

we see that ${\cal F}_{\ast{\bullet}}$ vanishes in the third order:

[TABLE]

Note also that Eqs. (116), (120) and (124) illustrate self-consistency check for the background-field condition (108).

One can continue and prove by induction that ${\cal F}_{\ast{\bullet}}$ vanishes in an arbitrary order in ${\bar{G}}_{\ast{\bullet}}^{n}$ and therefore the field ${\cal A}_{\mu}$ is a pure gauge

[TABLE]

Now we shall demonstrate that the matrix $\Omega$ satisfies our requirement (35). Since $C_{\ast}(x_{\ast}\rightarrow-\infty,x_{\bullet})~{}=~{}0$ due to Eq. (106), we get

[TABLE]

Similarly,

[TABLE]

One can also construct the expansion of matrix $\Omega$ in powers of ${\bar{A}}_{\ast}$ and ${\bar{B}}_{\bullet}$. For example, up to the fifth power of the ${\bar{A}}_{\mu}$ fields

[TABLE]

Now, for each $x_{\perp}$ we solve auxiliary 2-dimensional classical problem (105) and find

$\Omega(x_{\ast},x_{\bullet},x_{\perp})$ satisfying the requirement (35).