Improved evolution equations for QCD

Abstract Evolution equations based on probabilistic interpretation of QCD are derived. The main improvement is the inclusion of exact off mass shell kinematics in parton branching, which permits a realistic Monte Carlo simulation of the jet evolution.

In perturbative QCD the scaling violating Q2.depen_ dence of structure and fragmentation functions formally arises as a result of summation over all orders in leading logarithmic approximation (LLA) [1 ]. In axial gauge the leading contribution is obtained by restricting the sum over the tree diagrams with dressed propagators and vertices [2], while interference between different diagrams can be neglected. A consequence of these approximations is that the evolution of partons can be interpreted probabilistically [3], in fact the formation of a jet of partons initialized by a single parton resembles a shower due to an electron or a photon penetrating matter [4]. Konishi, Ukawa and Veneziano [5] developed a simple set of rules -Jet calculus -to derive expressions for inclusive multiparton distributions and parton multiplicities based on LLA. Recently Bassetto et al. [6] extended some of these considerations to study e.g. transverse momentum effects for multiplicity distributions as well as preconfinement [7].
The formulas derived on the basis of LLA are elegant and relatively simple. However, the simplicity is due to strong kinematic approximations and there are indications from, e.g. QED, that LLA can be quantitatively rather misleading. The main purpose of the present work is to derive rules and quantities which preserve exact kinematical relations in the parton branching process. In this way one can estimate the size of cor-rections to LLA. However, we have not really improved the QCD treatment beyond LLA. Our method is also suited to an efficient numerical calculation of any exclusive or inclusive quantity at the parton level using Monte Carlo simulation.
The kinematics relations are simplest when one uses as the scaled variable the light cone variable z denoting the energy + longitudinal momentum fraction of a parton from the parent parton: In the branching vertex p -+ p 1 + P2 the following exact formula holds where PT is the transverse component ofp 1 perpendicular to p. In LLA the dominant contribution comes from p2, _2 p2 the region 1 ta2 '~ , thus p2 and z-distributions become independent. In the phase space integral pl 2 and p2 are integrated up to p2, which results in a numerically significant contribution from a kinematically forbidden part of the phase space.
Let us next consider the discontinuity of a quark propagator or, in probabilistic interpretation the distribution of the valence quark of (off shell) mass Qo in the fragments of a quark created e.g. in e+e -annihilation at c.m. energy Q >> Qo" As in LLA [1 -6] we shall consider tree diagrams only, further we shall replace the dressed vertices and propagators by the running coupling constants and bare propagators, while we keep the exact kinematics relations and the correct phase space integral. Remembering that the probability for finding a valence quark in a quark is unity, we obtain Pm(Q2;q21 .... q22,p2 ) where P0 = 1 and (see fig. 1 (s) For clarity we have assumed that gluons do not decay and in order to get a finite result we take q2 = Qg. The factor Nq (Q2) is the unbalanced renormalization constant which is left over, when dressed propagators aqd vertices are recombined into running coupling constants [2]. Alternatively we may consider eq. (3) as a definition of Nq(Q 2) according to the Kinoshita-Lee-Nauenberg theorem. Physically interpreted Nq(O 2) is the probability that the quark is created (in the e+eannihilation vertex or some other hard process) directly with mass Qo and hence does not emit gluons corresponding to the term U0 in eq. (3). To generalize eq. (3) to contain gluon decays we have to sum over all possible branchings of each of the gluons. If there were no kinematic correlations between gluons and the valence quark under consideration in eq. (3), each branch would just contribute a factor 1. In the absence of such correlations there exists thus a renormalization condition for the gluon "form factor" Ng(Q 2) similar to eq. (3). Nq and Ng are the QCD analogues of the Sudakov form factor. Due to the cutoff O0 they are both finite, but will rapidly tend to zero, when Qo decreases. In the case of QCD eq. contains c~ s (4Q2), and the condition c~ s (402)/7r < 1 must be satisfied, hence O0 > A/2. The normalization condition (3) and the limitation to the class of tree diagrams specify the solution completely without further approximations. It is simplest to represent the solution in the form of an integral equation as follows: Firstly we shall redistribute the factors Nq and Ng in such a way that the propagator for a parton of type a becomes Na(p2 ) 2 2 O~s(Pa)dP a where p2 is the off shell mass squared of the parent patton. As in eq. (3) at each vertex a ~ bc there will be a factor where Pa~bc are suitably normalized Altarelli-Parisi decay functions without 6-functions. Sum of the tree diagrams using rules (6 -7) gives directly the cross section. Without loss of generality we can assume Nq(Q 2) = Ng(O 2) = 1, the normalization condition (3) can then be written in the form Q2 Na(O2) 1 =Sa(Q2 ) + f Na(p2--~ Pa(p2)dp 2 ,

Nc(p )
as the final result. The indicated sum applies only for gluons, in that case Pg = Pg~gg + Pg_.q~ corresponding to the two decay channels. We emphasize that the equation obtained does not involve any approximations after the form of the matrix elements in terms of the Altarelli-Parisi functions and renormalization conditions (8) have been fixed. For a classical relativistic branching process this would be an exact solution.
To assess the significance of eq. (9) it is important to see how it reduces to the known results in LLA. In view of the fact that we have not really rigorously solved QCD beyond LLA but only indicated a possible form of the generalization based on branching ideas, the fact that eqs.

Xn( @ o0 2-5
If one could set a s = (~s (Q0 2) this would agree with the Sudakov form factor derived in ref. [8]. In our case, however, the steps leading to eq. (1 l) cannot be numerically justified. Among others, notice that eq. (11) does not satisfy the factorization property discussed below.
In order to carry out a numerical solution of eqs. (8)-(9) we first note that eq. (8) simply tells that F(Q 2) is the logarithmic derivative of N(Q 2) Fa(Q2 ) = -dO ~ log Na(Q2 ) QZ Na(O2)=exp{_ f pa(p2)dp2} " (12) Substituting this into eq. (9) one obtains coupled nonlinear integral equations for I'a(p2), which are rather easy and fast to solve numerically. The simplicity is due to the fact that Pa(p 2) is given in terms of pa,(p~), where x/~ ~< @-7 _ Q0, hence the solution can be integrated in a single iteration. Values of Nq(O 2) obtained from this procedure are plotted in fig. 2. One can see that a change of Qo has a noticeable effect only at small virtual masses, which indicates the correct infrared behaviours. The shape of Ng(O 2) is very similar, except that the slope is approximately a factor 2 steeper than that of Nq(Q2). Fig. 2 also shows for comparison the LLA result based on eq. (1 l) [8]. Tile right hand side of eq. (11) describes a transition from the mass scale Q2 to Q2, which should be given by the ratio of N-functions. For comparison (in accordance with eq. (22) of ref. As the kinematical constraints are correctly included it is straightforward to proceed to Monte Carlo simulation of parton branching. One simply has to generate three numbers z, p2, p2 at a given branching vertex according to the independent distributions P(z), separately Pg.ogg and lPg~q~ as well as Ng~gg, Ng~q{ defined in an obvious way using eq. (12), since the generation of the gluon mass determines at the same time the future decay channel. The above Monte Carlo simulation nicely demonstrates the probabilistic interpretation of our formalism. Once a parton has been created with definite momentum (and for gluons with the information about the decay channel) its decay is independent of what happens to the other partons. In LLA the distributions of z, p2 and pc 2 in a vertex a -+ b,c would also be independent, which is impossible due to kinematics. In our case these distributions are otherwise independent except that they are coupled by the 0-function. The normalization conditions (8), which are crucial for the correct cancellation of infrared real and virtual singularities, together with the integral equations (9) take care of the proper normalization of the joint probability density. A quantity like Nq(O 2) is interpreted as the probability for a quark created fiom a parton of mass O to get the mass Oo without kinematics constraint. Notice that with probabilistic interpretation all parton form factors describing a transition (6) Q2 _+ pa 2 will be factorizable Na(Q2)/Na(p2a), since the probability N(O 2 -+ 002) (=-Na(O2)) = N(O 2 _+ p2) X N(p 2 -+ 002). Correspondingly the total gluon form factor Ng(O 2) is the product Ng = Ng_.+gg X Ng+q~, since the decay rates are additive Pg = I~g_~gg + Pg._.q~. In order to generate a gluon mass in practice one has to generate two masses using Ng_+gg and Ng__.q~ independently and choose the larger mass. This is of course distributed according to [8 (p2 _ O2) + Fg(p2)]/Ng(p2). In fig. 3 with Q02 = 6 GeV 2. It is evident that the effect of exact kinematics is to slow down the evolution appreciably, hence the standard estimates based on LLA can at best be qualitatively correct. Notice, however, that in our case the evolution ends at 4002(4 GeV 2 in fig. 3), since Na(Q2 ) = 1 for Q2 ~< 4Qo 2. in fig. 4 we have made a comparison with data on opposite side correlation from PETRA as measured by PLUTO [10] at 30 GeV. For this purpose we have included fragmentation into hadrons using a string picture [11 ]. Firstly we arrange all partons created in the evolution into preconfinement clusters [7]. After that each gluon is further split into qg:l using Pg_+q~ distribution: The resulting strings are split into hadrons with a Feynman-Field type algorithm [12] keeping the exact energy-momentum conservation for each string. Opposite side correlations should be sensitive to small deflection of quark and the antiquark due to nearly collinear gluon emissions. As pointed out by the PLUTO collaboration [10] as standard Monte Carlo program without QCD evolution was in slight disagreement with data on angular correlations. With our present formalism the QCD evolution seems to be  sufficient to produce the required effect in rough agreement with data.
In conclusion, we have presented a simple way of estimating the subleading log corrections to QCD evolution equations and shown that these corrections are likely to be sizeable. The method of generalizing QCD evolution equations leads in a straightforward way to Monte Carlo simulation of parton branching. First applications indicate that agreement with data is improved; the method can be applied to several other processes.