Summations of large logarithms by parton showers

We propose a method to examine how a parton shower sums large logarithms. In this method, one works with an appropriate integral transform of the distribution for the observable of interest. Then, one reformulates the parton shower so as to obtain the transformed distribution as an exponential for which one can compute the terms in the perturbative expansion of the exponent. We apply this general program to the thrust distribution in electron-positron annihilation, using several shower algorithms. Of the approaches that we use, the most generally applicable is to compute some of the perturbative coefficients in the exponent by numerical integration and to test whether they are consistent with next-to-leading-log summation of the thrust logarithms.

The logarithms L j (r) arise in QCD from the soft and collinear singularities of the theory. These same soft and collinear singularities are contained in the splitting functions of a parton shower algorithm. Thus running a parton shower event generator to calculate σ J (r) will produce an approximation to the series in Eq. (1). That is, the parton shower approximately sums the large logarithms. The object of this paper is to investigate the form of the result of this summation. 1 To exhibit the summation of logarithms, we rearrange the parton shower algorithm so that it is specialized to calculate just σ J (r) and so that it expresses σ J (r) directly in terms of an exponential The integral of S Y (µ 2 ; r) in the exponent has an expansion 2n j=0 e(n, j) L j (r) . (3) The operator S Y (µ 2 ; r) is determined by the parton splitting operator S(µ 2 ) in the original shower. This gives one direct access to the coefficients e(n, j). With this representation, one has the potential to prove that e(n, j) = 0 for j > n + 1. The terms with j = n + 1 are called leading-log (LL) terms and the terms with j = n are called next-to-leading-log (NLL) terms. One also has the potential to prove that e(n, j) for j = n + 1 and for j = n are what is expected in full QCD if a full QCD result is known. Our plan in this paper is to develop the general theory behind the representation (2) along the lines of Ref. [1]. In this exposition, we also present the main steps of the construction of Ref. [1] in a form that, in our opinion, makes these steps more transparent. In a companion paper [2], we apply the representation (2) to an important example, the thrust distribution in electron-positron annihilation. We consider just the thrust distribution and not other distributions involving large logarithms. However, we look in some detail at how the exact form of the shower algorithm affects the results.

II. PARTON SHOWER FROM PERTURBATION THEORY
The starting point is the perturbative cross section for an infrared safe observable in hadron-hadron collisions. We describe this briefly here. A more detailed explanation can be found in Ref. [1].
The parton shower is described using operators on a vector space, the "statistical space," that describes the momenta, flavors, colors, and spins for all of the partons created in a shower as the shower develops. The colors and spins are quantum variables and are described using a density matrix. We use this description in the parton shower event generator Deductor [3][4][5][6][7][8][9][10]. The general theory includes parton spins, so we include spins here even though Deductor simply averages over spins. With m final state partons plus two initial state partons with labels "a" and "b," the partons carry labels a, b, 1, 2, . . . , m. The partons have momenta {p} m = {p a , p b , p 1 , . . . , p m } and flavors {f } m . We take the partons to be massless: p 2 i = 0. For color, there are ket color basis states |{c} m and bra color basis states {c ′ } m |. We use the trace basis, as described in Ref. [3]. For spin, there are ket basis states |{s} m and bra basis states {s ′ } m |. Then the m-parton basis states for the statistical space are denoted by |{p, f, c, c ′ , s, s ′ } m ). A vector |ρ) in the statistical space is a linear combination of the basis states.

A. Perturbative cross section
If the QCD matrix element is calculated up to a given order, α K s , the cross section is Here the renormalized perturbative QCD density operator is represented by a vector in the statistical space |ρ(µ 2 r )). It is based on the exact matrix elements and contains all the possible partonic final states at order K. The density operator is already renormalized, typically in the MS scheme, thus it is independent of the renormalization scale, µ 2 r , up to the desired order The next factor in Eq. (4) is the operator of the bare parton distribution functions (PDFs), Here the circles, a • b, represent convolutions in the momentum fraction variables. The renormalized PDF operator for the hadron-hadron initial state is F r (µ 2 r ). The corresponding MS subtraction of initial state singularities is done by the Z F (µ 2 r ) operator, which contains factors 1/ǫ n in dimensional regularization. As described in Ref. [1], one should typically use something other than the MS scheme to define the parton distribution functions used internally in the shower. The factor K(µ 2 r ) transforms to the shower scheme for the parton distribution functions F r (µ 2 r ). The bare PDF is scale independent, This equation leads to the proper evolution equation of the renormalized PDFs. The next factor in Eq. (4) is the operator O J (r) representing an infrared (IR) safe measurement, characterized by a set of parameters r.
After applying these operators, we have a sum and integral over basis states |{p, f, c, c ′ , s, s ′ } m ). Finally, we multiply by the statistical bra vector (1| and obtain a cross section after performing the integrations using (8) (The spin states |{s} m are orthogonal and normalized, but the color states |{c} m in the trace basis that we use are not orthogonal and some of them are not normalized exactly to 1 [3].) If the calculation includes perturbative contributions up to α K s , then there is an error term O(α K+1 s ) in Eq. (4). The formula is based on standard QCD factorization for infrared safe observables. This has power suppressed corrections of order Λ 2 QCD /Q 2 J (r) where Q 2 J (r) is the lowest scale that the measurement operator O J (r) can resolve. In the rest of this paper, we mostly omit explicit mention of these error terms.
The expression in Eq. (4) simplifies substantially in electron-positron annihilation. In this case, we can replace the operator F 0 by 1.
We point out that Eq. (4) is valid only in d = 4 − 2ǫ dimensions. It is not directly useful for practical calculations.

B. IR singular operator
To define a good subtraction scheme for a fixed order calculation one can use the IR singular operator D(µ 2 r , µ 2 s ) [1]. This operator has a perturbative expansion The operators D (n) (µ 2 r , µ 2 s ) are key to defining a parton shower algorithm in a general framework. For a first order shower, one uses only D (1) (µ 2 r , µ 2 s ), but in a general framework we consider D (n) (µ 2 r , µ 2 s ) for any n. This operator describes the IR singularity structure of partonic states |ρ(µ 2 r )). When D (n) (µ 2 r , µ 2 s ) acts on a state {p, f, c, c ′ , s, s ′ } m it produces new states {p,f ,ĉ,ĉ ′ ,ŝ,ŝ ′ }m with m ≤m ≤ m + n such that the IR singularities of match the singularities of nth order QCD Feynman diagrams that connect these two states. Here the singularities include the factors 1/ǫ from virtual loop diagrams and they include the singular behavior of the diagrams when any two or more momentap become collinear or some of thep become soft. The operator D (n) (µ 2 r , µ 2 s ) depends on two scales, the standard renormalization scale µ 2 r and the shower scale µ 2 s . The shower scale acts as an ultraviolet (UV) cutoff that separates the IR and UV regions. All IR singularities are included, but only regions near the singularities with scales k 2 satisfying k 2 < µ 2 s are included. There is, of course, some freedom in choosing how the division between IR and UV is defined. Different prescriptions lead to differences in the shower ordering prescription in the parton shower algorithm produced by D (n) (µ 2 r , µ 2 s ). The singular operator is based on the MS renormalized matrix elements and is independent of the renormalization scale. Thus we have This allows us to choose the renormalization scale conveniently.
In order to avoid large logarithms of µ 2 r /µ 2 s , it is useful to relate the renormalization scale to the shower scale. We define Then we can avoid large log(κ r ) factors by choosing κ r of order 1. The singular operator is perturbative and we can always define its perturbative inverse operator, by working order by order in the perturbative expansion.

C. Fixed order cross section
We can make Eq. (4) more useful by inserting 1 in the form DD −1 , We notice that the expression D −1 (µ 2 r , µ 2 s )|ρ(µ 2 r )) is well defined in d = 4 dimensions since the inverse of the singular operator removes all the IR poles of |ρ(µ 2 r )). Accordingly, we define the subtracted hard matrix element by This gives us We will use Eq. (15) to explore parton showers. First, however, suppose that we are interested only in the fixed order cross section. Then we can choose the scale µ 2 s small enough that the measurement operator O J (r) does not resolve parton momentum scales of order µ 2 s . Then One can calculate (1|F 0 D(µ 2 r , µ 2 s ) in d = 4 − 2ǫ dimensions. The operator D(µ 2 r , µ 2 s ) creates singularities, but the initial state singularities are removed by the operator Z F (µ 2 r ) in F 0 and the final state singularities cancel after we multiply by (1| and integrate over the parton variables. Thus we obtain a finite result in the ǫ → 0 limit.

D. Operators V and X1
The operators D(µ 2 r , µ 2 s ) and F 0 are defined only in d = 4 − 2ǫ dimensions and are singular as ǫ → 0 and as parton momenta become soft or collinear. However, we have noted that (1|F 0 D(µ 2 r , µ 2 s ) is finite in d = 4 dimensions. It will prove useful to introduce an operator, V(µ 2 r , µ 2 s ), that is finite in four dimensions, does not change the number of partons, leaves the parton momenta and flavors {p, f } m unchanged, and satisfies The operator V(µ 2 r , µ 2 s ) leaves {p, f } m unchanged, but it can act non-trivially on the color and spin space. Eq. (17) does not fully define the color and spin content of V(µ 2 r , µ 2 s ). We discuss the definition further in Sec. IV, but for now, we need only Eq. (17).
Using V(µ 2 r , µ 2 s ) we define a singular operator X 1 (µ 2 r , µ 2 s ) as The "1" subscript distinguishes the operator X 1 from the operator X used in Ref. [1] and suggests the normalization condition (19).
With V and X 1 , the cross section in Eq. (15) can be written as This form will be useful to help us define a parton shower. Before we continue with the discussion of the parton shower cross section we introduce a more compact notation for operators with renormalization scale dependence. According to Eq. (11) the renormalization scale is always related to the shower scale, thus we can define The PDF operator depends only on the renormalization scale and in this case the convention is a little different, The functions specified above then depend on κ r , but we do not display this dependence. With this more compact notation, Eq. (20) is written as

E. Operator U and parton shower
The formula for the cross section σ J given in Eq. (23) is of limited usefulness if the scale Q 2 J (r), representing the lowest scale that the measurement operator O J (r) can resolve, is much smaller than the scale µ 2 h of the hardest momentum transfer in ρ h (µ 2 ) . When that happens, σ J will contain logarithms log(µ 2 h /Q 2 J (r)) that need to be summed by looking for the most important terms at all orders of perturbation theory. To that end, one can use a parton shower algorithm.
To provide a parton shower, first set the scale µ 2 in Eq. (23) to µ 2 h . Then define a scale µ 2 f that is certainly smaller than Q 2 J (r). Typically, one chooses µ 2 f on the order of 1 GeV 2 . Finally, Since µ 2 f < Q 2 J (r), the operator O J (r) does not resolve partons at the scale µ 2 f . Thus O J (r) commutes with X 1 (µ 2 f ), giving us With the use of Eq. (19), this is The operator X −1 This is the shower operator. It generates a parton shower starting at the scale µ 2 h and ending at the scale µ 2 f . Because of Eq. (19), the shower operator is probability preserving Using the notation U(µ 2 f , µ 2 h ), the cross section is (29) We have perturbatively calculated matrix elements with their IR divergences subtracted in ρ h (µ 2 h ) . Then the operator F (µ 2 h ) supplies parton distribution functions. The factor V(µ 2 h ) serves to sum threshold logarithms [1,11]. An approximation to this factor is contained in Deductor although it is lacking in other current parton shower event generators. Next, the operator U(µ 2 f , µ 2 h ) generates the parton shower and the operator O J (r) measures the desired observable in the multiparton state created by the shower. Finally, we multiply by (1| and integrate to get the desired cross section. We discuss U(µ 2 f , µ 2 h ) and V(µ 2 h ) in more detail in Secs. V and VI.

III. OBSERVABLE DEPENDENT SHOWER EVOLUTION
The operator O J (r) in Eq. (29) could represent any infrared safe observable. In this paper, we have a particular sort of operator in mind. Consider, for example, the transverse momentum distribution of a Z boson produced in the Drell-Yan process. The operator that measures the transverse momentum k ⊥ of the Z boson is defined aŝ where k Z ({p} m ) is the transverse momentum of the observed Z boson. The standard method for summing logarithms of k 2 ⊥ /M 2 Z is to start with the Fourier transform of the k ⊥ distribution. To measure this with a parton shower event generator, we can use the measurement operator We let O Z (b) serve as an example of the observable O J (r) that we consider in this paper. There are many other similar examples. We will need one property of the observable O J (r) beyond infrared safety: we assume that the operator O J (r) has an inverse O −1 J (r). To analyze the cross section σ J (r), we start with the representation (23) with µ 2 = µ 2 h , Define an operator Y(µ 2 ; r) that is finite in d = 4 dimensions, leaves the number of partons and their momenta and flavors unchanged, and is related to X 1 by Then define a new version of X 1 that depends on the measurement parameters r by This gives us and Then our cross section is With the use of Eq. (35), and commuting O J (r) past V(µ 2 h ) and F (µ 2 h ), which do not change the partonic state, this becomes Here we measure O J (r) at the hard state ρ h (µ 2 h ) , obtaining typically a very simple result. Then we measure O J (r) inside the operator Y(µ 2 h ; r). This operator has the potential to sum large logarithms.
We can also relate Y(µ 2 ; r) to the shower operator U(µ 2 f , µ 2 ) with a small final scale µ 2 f . From Eq. (33), we have Since µ 2 f < Q 2 J (r), the operator O J (r) does not resolve partons at the scale µ 2 f . Thus O J (r) commutes with X 1 (µ 2 f ), giving us Recall from Eq. (19) that (1|X 1 (µ 2 ) = (1|. This gives us That is, we compare two calculations. In the first calculation, we generate a parton shower down to a very small scale starting with any statistical state at a scale µ 2 . Then we measure O J (r) inclusively using (1|O J (r).
In the second calculation, we first operate with O J (r) on the state at scale µ 2 then measure Y(µ 2 ; r) inclusively using (1|Y(µ 2 ; r). These two calculations give the same result.

IV. THE OPERATOR MAPPING P
In Sec. II D we defined an operator V(µ 2 ) which is to obey Eq. (17), 1 V(µ 2 ) = 1 F 0 D(µ 2 )F −1 (µ 2 ). In Sec. III, we defined an operator Y(µ 2 ; r) in the same way. In each case, we start with a singular operator A and we want to define a second, nonsingular, operator B with the property When the operator B acts on an m-parton basis state {p, f, c, c ′ , s, s ′ } m , it is to leave the number of partons, their momenta, and their flavors unchanged. It may, however, act non-trivially on the colors and spins. These requirements do not fully specify B. We can be somewhat more definite by requiring that there be a linear mapping A → B, which we write in the form This mapping must satisfy and [A] P must leave m and {p, f } m unchanged, The requirement (43) is then a restriction on the spin and color matrix A, We can place another requirement on [· · ·] P : if A has the property that it leaves m and {p, f } m unchanged, then One consequence of this is that [[A] P ] P = [A] P . These requirements do not fully specify the mapping [· · ·] P . For this paper we do not need to be more specific. However in Ref. [2] we provide an example (without spin) that is useful for the analysis of a first order shower.
We will find that the combination A − [A] P appears frequently in formulas. It useful to define an operation [· · ·] 1−P by (49)

V. GENERATOR OF SHOWER
We now turn to a more detailed study of the operator U(µ 2 f , µ 2 h ) that creates a parton shower between a hard scale µ 2 h and a small, cutoff scale µ 2 f . The generator of this shower evolution is the operator Here, we differentiate with respect to the shower scale. Because of Eq. (10) (with the use of Eqs. (17) and (18)), this is the same as Because of Eq. (35), Eq. (51) gives us a differential equation for U We use the notation to represent the solution of this equation. Here T indicates the instruction to order the operators S(µ 2 ) with the smallest µ 2 to the left.

VI. THE THRESHOLD FACTOR
In Sec. II D we have defined an operator V(µ 2 r , µ 2 s ). With our notation in Eq. (21) for the scale dependence of V, the crucial property given in Eq. (17) can be written In Eq. (29) or Eq. (38), the perturbative expansion of V(µ 2 h ) contains large logarithms [1,11,12]. These are the much studied threshold logarithms [13]. We sum the threshold logarithms by writing V(µ 2 h ) as an exponential. Define Then V(µ 2 h ) can be written as Define a generator operator S V (µ 2 ) by Then U V (µ 2 2 , µ 2 1 ) is the solution of the differential equation We write the solution of this equation as As long as we expand the running coupling α s in Eq. (60) to some finite order in α s (µ 2 h ), the integral in Eq. (60) is convergent in the limit µ 2 f → 0 [1]. Thus V(µ 2 ) at small scales is almost the unit operator, That is

VII. PERTURBATIVE EXPANSIONS
The operator S(µ 2 ) can be expanded in powers of α s (µ 2 r ) = α s (κ r µ 2 ): In the general theory from Ref. [1], S(µ 2 ) is constructed from the singular operator D(µ 2 r , µ 2 s ). If we use only the first order part D (1) (µ 2 r , µ 2 s ) of D because that is all we know, then all we get is S (1) (µ 2 ). However, in a practical parton shower program (such as the Λ-ordered Deductor), one often takes a guess at approximate higher order contributions S (n) . The approximate form is obtained by changing the argument of α s in the splitting functions to κ r k 2 T and, additionally, making a special choice for κ r . Expanding α s (κ r k 2 T ) in powers of α s (κ r µ 2 ) then produces contributions S (n) (µ 2 ) for n > 1.
In Deductor, the first order contribution has three parts [1,11]: The operator S (1,0) (µ 2 ) describes parton splitting, changing an m parton state to an m + 1 parton state. The operator [S (1,0) (µ 2 )] P leaves m and {p, f } m in an m parton state unchanged, although it can modify the color state. 2 In a leading color parton shower, the color is unchanged and the eigenvalue of this operator then gives the order α s contribution to the integrand in the exponent of the Sudakov factor that represents the probability not to split between two scales. The final operator, S (0,1) iπ (µ 2 ), leaves m and {p, f } m unchanged. It gives the imaginary part of virtual graphs [1,11] and obeys (1|S The operator S V (µ 2 ) has a perturbative expansion The first order operator S V (µ 2 ) has the form [1,11,12] S (1) Here [S (1,0) (µ 2 )] P is proportional to the integral of the first order splitting function over the splitting variables and appears also in Eq. (64). In the third term, [F (µ 2 ) • P (1) ] denotes the convolution of F (µ 2 ) with the first order PDF evolution kernel P (1) . In the second term, is the derivative with respect to the shower scale of the singular operator for a one loop virtual graph. It is sometimes assumed that the effect of virtual graphs and PDF evolution cancels the integral over the splitting variables of parton splitting [14]. However, this cancellation is not complete, so that the effect of S V (µ 2 ) is quite important [12,14].

VIII. GENERATOR OF Y
We now turn to a more detailed study of the operator Y(µ 2 ; r). This operator sums logarithms, so we want to write it as an exponential. Define This gives us a differential equation for Y(µ 2 ; r) Using Eq. (75) then gives us The operators Y and S Y are nonsingular operators that leave the number of partons and their momenta and flavors unchanged. Thus we can use the mapping [· · ·] P defined in Sec. IV to write this as The expansion of Y in powers of α s starts at Y = 1 + · · · , so a useful way to write this is Now we can use Eqs. (80) and (73) recursively to generate S Y and Y in powers of α s . We write with For S Y , Eq. (80) gives This gives us S (n) Y for j < n and Y (k) for k < n.
For Y we use Eq. (73), in which an integration over an intermediate scaleμ 2 appears. We can expand α s (κ rμ 2 ) in powers of α s (κ r µ 2 ) in the form with coefficients γ derived from the QCD β-function. Using this expansion in Eq. (73), we obtain This gives us Y (n) if we know Y (k) for k < n and S (j) Y for j ≤ n.
These recursion relations successively generate S Y , Y (2) , . . . . The first order terms are and We now outline how the operator Y(µ 2 ; r) can be used. This operator is the key to calculating an observable cross section σ J (r) according to a parton shower algorithm. The operator O J (r) that defines this cross section must be infrared safe. That is, there is a scale Q 2 J (r) such that σ J (r) does not resolve parton splittings at scales µ 2 smaller than Q 2 J (r). In order to define Y(µ 2 ; r), the inverse operator O −1 J (r) must exist. The anticipated use case is that there is a distribution of direct interest that involves large logarithms and the logarithms can be summed analytically by taking an integral transform of the distribution that depends on parameters r. Then σ J (r) represents the value of this integral transform. In a companion paper [2], we examine an important example, the thrust distribution in electron-positron annihilation. Then one uses the Laplace transform of the thrust distribution and r is the Laplace parameter ν.
In the applications that we have in mind, the perturbative expansion of σ J (r) contains powers of a large logarithm L(r) when the parameter or parameters r approach some limit. Typically, we have In favorable cases, there is an analytical formula that sums these logarithms in the form It is crucial here that the maximum power of L at order α n s is j = n + 1, not 2n. We can say that a σ J (r) with this property exponentiates. One never knows all of the coefficients d(n, j), but when the coefficients for j = n + 1 are known, we can say that the formula sums the logarithms at the leading-log (LL) level. When the coefficients for j = n are also known, we can say that the formula sums the logarithms at the next-to-leaadng-log (NLL) level.
In some important cases, the color space for the partons involved in the hard scattering process is trivial. For instance, for shape observables in electron-positron annihilation, there is only one color basis vector for the qq state in e + e − → qq. Then the coefficients d(n, j) are numbers. The initial partonic state in hadron-hadron scattering has a nontrivial color structure. Then the coefficients d(n, j) may be integrals of matrices in the parton color space, with some specification for the ordering of noncommuting matrices in the exponent.
What does a parton shower algorithm say about σ J (r)? Different parton showers can give different answers, so we should have a particular parton shower algorithm in mind.
We have seen that there are two ways to express σ J (r) as given by a parton shower. First, we can use Eq. (29), (90) Typically the splitting operator S in U(µ 2 f , µ 2 h ) is based on lowest order perturbation theory, as discussed at the beginning of Sec. VII. Additionally, V(µ 2 h ) is present in Deductor, but for many parton shower algorithms V = 1. Equation (90) says to run the parton shower to its cutoff scale and then measure the observable by applying (1|O J (r). The perturbative expansion of this result has the form (88), but not directly the form (89). One can run the corresponding parton shower event generator to obtain a numerical result with statistical errors and other numerical errors. Even with errors, it is possible [2,15] to use numerical results from Eq. (90) to check these results against a known QCD analytic result.
The second way to express σ J (r) as given by a parton shower is contained in Eq. (38), The operator S Y is obtained from the shower generator S using Eqs. (83) and (85). This second expression for σ J (r) gives exactly the same σ J (r) as given by Eq. (90). However, now the logarithms L appear in the exponent in S Y . Thus we have a representation that is very close to the representation in Eq. (89). The exponent in Y is 3 If we use the perturbative expansion of S Y , this is I(r) = as numerical integrals and check how many powers of L(r) they contain. Although one can never check every term in ∆I(r), this method has the advantage that if the check for NLL summation fails for any one contribution, then we know that NLL summation fails. In Ref. [2], we present the analysis of I(r) in a form that is somewhat less general than the analysis of this paper but is better adapted to practical applications. Then we analyze I(r) analytically and numerically for the trust distribution in electron-positron annihilation.

X. CONCLUSIONS
It is, we think, of some importance to understand how accurately a parton shower algorithm sums large logarithms in an observableσ J (v).
In analytical approaches to summing such logarithms, one typically defines an integral transform of the original distribution so that one considers a cross section σ J (r) that depends on parameters r. Then the perturbative expansion of σ J (r) contains large logarithms L(r).
Sometimes, one can compare the results of the shower for σ J (r) to the results in full QCD by writing the same differential equations as for full QCD but applying the differential operators to the shower approximation rather than full QCD [16,17]. This method has he disadvantage that one needs a separate and quite elaborate analysis for each observable to be studied.
An alternative is to calculate the observableσ J (v) numerically with the parton shower event generator of interest and to compare the result with a known QCD result [2,15]. This method can work, at least for electron positron annihilation, but presents significant numerical challenges.
We have presented a reformulation of the calculation of σ J (r) according to a parton shower so that the large logarithms appear directly as an exponential. The exponent can be expanded perturbatively. This gives us a path to an analytical understanding the summation of these logarithms in the parton shower. It also provides a simple way to test this summation numerically. In a companion paper, we find interesting results [2] for the thrust distribution in electron-positron annihilation.

ACKNOWLEDGMENTS
This work was supported in part by the United States Department of Energy under grant DE-SC0011640.