NLLjet: a Monte Carlo code for e+e- QCD jets including next-to-leading order terms

Computer Physics Communications 64 (1991) 67—97 North-Holland

67

NLLjet: a Monte Carlo code for e~e QCD jets including next-to-leading order terms Kiyoshi Kato Kogakuin University, Department of Physics, Nishi-Shinjuku 1-24, Shinjuku, Tokyo 160, Japan

and Tomo Munehisa Yamanashi University, Faculty of Engineering; Takeda 4-3, Kofu 400, Japan Received 20 March 1990; in revised form 8 August 1990

NLLjet, a Monte Carlo generator for the parton shower in eke-annihilation has been developed. It can be used to determine the QCD A-parameter from eke-jet events since it includes the next-to-leading order contribution.

PROGRAM SUMMARY Title of program: NLLjet

No. of lines in combined program and test deck: 8572

Catalogue number: ABZA

Keywords: perturbative QCD, jets in e~e-anthhilation,parton shower, Monte Carlo simulation

Program obtainable from: CPC Program Library, Queen’s University of Belfast, N. Ireland (see application form in this issue) Computer for which the program is designed and others on which it is operable: FACOM, HITAC, IBM and others with a FORTRAN77 compiler Computer: FACOM-M780; Installation: National Laboratory for High Energy Physics (KEK), Tsukuba Ibaraki 305, Japan Operating system: OSIV/F4 MSP

Nature of physicalproblem For the study of jets in high-energy reactions, a Monte Carlo event generator of QCD is indispensable. In order to determine the QCD fundamental parameter, i.e. the strong coupling constant, the parton-shower method including next-to-leading order terms is required. Method of solution Based on the formulation of the QCD parton shower beyond the leading order, all terms up to the next-to-leading order are implemented into NLLjet. Detailed discussion on the kinematics of the parton shower is presented in this paper.

Programming language used: FORTRAN77 High speed storage required: 187 Kwords No. of bits in a word: 32

Typical running time On FACOM-M780, NLLjet generates 120 events/s of parton 2. When it is linked shower W= 91 GeV and Q~=1 GeVwith a set of typical with thefor LUND6.3 hadronization model parameters, 24 hadronic events are generated per second.

Peripherals used: Terminal or card reader for input, magnetic disk and line printer for output

0010-4655/91/$03.50 © 1991

—

Elsevier Science Publishers B.V. (North-Holland)

68

K. Kato, I Munehisa

/ NLLjeI: a

MC code fore ~e

QCDjets

LONG WRITE-UP 1. Introduction Among the studies of hadron jets, the most important subject is to check whether they are precisely described by QCD and, after its validity has been confirmed, to determine the fundamental parameters of QCD. The quark mass and the parameter for QCD-vacuum (9-vacuum) do not play any important role in high-energy perturbative QCD [1]. There is only one remaining parameter: the coupling constant of QCD which is denoted as a8 = g2/4’rr. Through the so-called dimensional transmutation it is equivalent to the QCD A parameter. The standard approach for simulating hadronic events in the e ±e -annihilation is divided into two steps. Firstly we generate a set of partons (quarks and gluons) according to the perturbative QCD calculation and secondly we convert them into observable particles using a hadronization model. Although the whole process should, of course, be described by the QCD Lagrangian, the short of knowledge regarding the solution of QCD in the soft sector forces us to employ a two-step treatment. These steps are separated by introducing an energy scale which is denoted as Q 0 hereafter. In the first step we need some approximations since we cannot sum all the diagrams of the perturbative QCD. There are two methods to do this: the matrix element method (ME) and the parton shower method (PS) [2]. In ME, we calculate the matrix elements by means of conventional perturbation. Because of the increasing technical difficulty, the calculation has been limited only up to order a~(and tree calculation of order a~). In PS, using the method based on the renormalization-group equation, we sum up the contributions of dominant terms from all order diagrams. The PS that keeps only the leading terms is called PS-LL (leading logarithmic approximation) and the one that keeps up to the next-to-leading terms is called PS-NLL (next-to-leading logarithmic approximation). In order to determine the QCD A, the model must satisfy two requirements. The first is that the model be constructed on the basis of the approximation that is capable of determining A. The second is that the effect of the hadronization be minimized. To remove this effect, we should take Q0 as small as possible. Q0 is about 1 GeV because the onset of non-perturbative effect appears around 1 GeV. As for the first requirement, PS-LL is not appropriate for the determination of QCD A since in LL one cannot fix the renormalization scheme [1,3]. The second requirement suggests that ME is not suitable since Q0, which corresponds to the minimum invariant mass of a pair of partons in ME, must be much larger than 1 GeV due to the existence of the mass singularity. Around 1 GeV, the value of the strong coupling constant calculated in LL differs significantly from that in NLL. This fact is also a shortcoming of PS-LL since it shows that the LL approximation is not sufficient around 1 GeV. Summarizing, among methods at hand, only PS-NLL has the potential for determining the QCD A. The development of the QCD parton shower has a long history. The first important work was the invention of the “jet calculus” by Konishi, Ukawa and Veneziano in 1979 [4]. Soon afterwards, Odorico and others formulated the method of PS-LL and also developed programs for generating a parton shower in the e~e-anthhilation and in other processes [5]. In 1984, following the results of the study on the soft gluon, Marchesini and Webber introduced the concept of angular ordering into PS-LL [6]. It gave an answer how to separate the overlapping of phase space of evolving partons. Taking these developments into account, Sjöstrand and Bengtsson released the program LUND version 6 for e~e-anthhilation which incorporated ME and PS with their string fragmentation [7]. In the LUND program, they also introduced into PS-LL the matrix element for e~e q~gto reproduce the hard three-jet events. The introduction of the three-jet cross section was unavoidable since the LL approximation is not correct in the large transverse momentum region. However, the theoretical base for this introduction was unclear because the three-jet cross section is only a portion of contributions in the NLL order. —‘

K Kato, T. Munehisa

/

NLL1eI: a MC code fore ~e — QCDjets

69

The theoretical framework for going beyond the LL was built by Konishi, Kalinowski, Scharbach and Taylor in 1981 [8]. They also calculated quantities in NLL for the quark sector. The calculation for the gluon sector was done by Gunion, Kalinowski and Szymanowski in 1984 [9]. Hence all parts have been ready for the construction of a Monte Carlo program based on the PS-NLL method. In 1987, the authors of the present paper showed how to construct a Monte Carlo model for PS-NLL, which is not just a trivial assembly of existing parts, but needs a detailed consideration of the kinematics [10,11]. In the conventional perturbative QCD analysis we treat the inclusive quantities, i.e. unobserved partons are integrated out. On the other hand, in the Monte Carlo approach we deal with the exclusive distribution, i.e. we must determine 4-momenta of all partons. We have developed a Monte Carlo program for eke-jets, called NLLjet. It should be noted that this program generates only the 4-momenta of partons, but it does not include the mechanism for hadronization. As already remarked, we expect that the model-dependence of hadronization is minimized. Consulting the appendices, users can easily generate e~e hadrons, connecting NLLjet with a hadronization model. The purpose of this paper is to show how the NLLjet can be formulated. Then we concentrate on the technical aspects and the kinematical details. Practical information for the users is also provided. In the next section, the formulation of NLLjet is briefly reviewed. Sections 3 to 11 give a description of the parts of NLLjet presented in section 2. In appendix A, we summarize basic notations used throughout this paper. In appendix B, an explanation on program modules of NLLjet is given. Description of common areas is shown in appendix C. The control parameters are explained in appendix D. In appendix E, we show some examples. —

2. Fonnulation of NLLjet The formulation of the parton shower beyond the LL is presented in ref. [8]. Using the off-shell unitarity equation of Polyakov for e e -annthilation, these authors write down all the matrix elements for e e partons diagramatically. The next step is to impose a projection operator on each propagator to separate the terms with mass-singularity from others which we shall call finite terms. Here the choice of gauge is essential. In the physical gauge, a pair of the same propagator gives rise to the mass-singularity. A diagram that contains an interference-type subdiagram gives less logarithms from the mass-singularity and such a non-leading part can be merged into an effective vertex. After repeated application of the procedure, the resulting diagrams have the following property: The effective vertices are finite and the diagrams are diagonal, allowing us to simulate the parton shower as a stochastic branching process. The basic ingredients of the PS-NLL are the following four items: +

+

- —*

(1) Effective propagator in the NLL. This corresponds to the non-branching probability (non-emission probability) TINE, which is calculated by integrating the Altarelli—Parisi function in NLL order. The non-branching probabilities for quark and gluon, i.e. the Sudakov factors [12], are given by j~dx a~(k2)Pqq(X,

1ThE(Q,

Q~)

H~E(Q,

Qfl =exp[_fQ2~f*dx~~

=

ex~[_

JQ2

K2)]~

(~P~(x~ K2)

+

~Pgq(x,

K2))];

the kinematical boundaries x ± are given in section 3. The Altarelli—Parisi function,

Pab(X),

(2.1) gives the

70

K. Kato, T. Munehisa

inclusive probability for a 1~ Pqq(X)=CF 1

X

b

—*

+

/

NLLjet: a MC code for e + e - QCDjets

anything. The explicit form in LL is [13]

Pgq(x)=T~[x2±(i_x)2],

,

(2.2)

Pgg(x)=2CA[l~X+ 1ff~+x(1_x)j. 2.

Those NLL order, which are found in ref. [10], are more complex on K (a q in or g) gives the probability that a parton a of virtuality Q2 turnsandto depend one of virtuality branching. In our formulation, K2 is P21 [14]. =

=

Q~) Q~without

H~E(Q2, =

(2) Two-body and three-body effective vertices. In NLL, vertex functions for two-body branching, V12~(x 1,x2) (x1 + x2 1), include higher-order terms. The following channels exist: q—*q+g, g-~g+g, g—’q+4. =

3~(x

Vertex functions for three-body branching, V~ which are absent in the LL, ie. q~q+g+g, g-*g+g+g,

1, x2, x1) (x1 + x2 + x3

=

1), generate new processes

q~q+q’+4’, g-*g+q+~j.

To denote the branchings p(x 1) a(x1) + b(x2) and p(x 1) a(x1) + b(x2) + c(x3), we write f~~(x1,x2) and V~(x1, x2, x3) since the parton p is unique for a given set of daughter partons. The last argument is sometimes omitted in the following. =

—~

=

—*

(3) Effective vertices that couple to the electro-weak current. In the present order, q4 and q4g are the primary objects produced by the annihilation of eke. In the lowest order, the q4g cross section is dø dz~dz2

a8 =

~~CF

___________

(1— z1)(1

—

z2)’

(2.3)

where ø~ is the lowest-order cross section for e e calculated by Gottschalk and Schatz [15]. +

— —~

q4. We have also included the order a~cross section

(4) Kinematics. This includes the angular ordering [6,14], the double-cascade scheme [11] and separation of jets, choice of variables and so on. 2) in the infinite momentum The(IMF). evolution thebranching, parton shower recorded by a pair variables, (x, Q is determined by x and Q2 frame At of each x is is always conserved. Theofmagnitude of PT through p-momentum conservation. The rotational angle of PT around the z-axis is denoted by In NLLjet, there is no dynamical condition for 4, and it is determined at random after the evolution has been terminated. In the following, we sometimes denote the momentum of a parton as a triplet of variables, ~.

(x,

Q2~PT)•

First we generate partons at the primary vertex which couples to electro-weak current in e e -annihilation. We have two types of primary vertex, q~and q4g. Either of these channels is selected at random as is described in section 5. In the q4 primary vertex, both quark and antiquark get x 1 and PT 0. In the q~gprimary vertex, the light-cone fractions of three partons, x,, and the transverse momentum of the quark, PT’ are determined according to a distribution. The transverse momentum of the antiquark is chosen to be 0 and that of the gluon is ~PT• Virtuality of the parton is not yet determined here. Details regarding these points will be found in sections 5, 6 and 7. +

=

=

K Kato, T. Munehisa / NLLJet: a MC code fore ~e — QCD jets

71

After generation of the primary vertex, we develop the parton shower iteratively. (a) There is a parton that has a probability to branch. It has a definite light-cone fraction but its virtuality is not yet determined. The virtuality of the parton is less than Q~by a kinematical limitation. (b) Using the Sudakov factor, we determine its virtuality Q2 in the range of Q~iax> Q2> Q~.If Q2 < Q~, it cannot branch, so that it is marked inactive with the virtuality Q2 Q~.If the determined virtuality is greater than Q~,it can branch and proceeds to the next step. (c) Using the Altarelli—Parisi function P(x), we determine the light-cone fraction, x ( x < x
—

~1~dx~2~

~

~Pgq(x)

and determine x. (d) In steps (b) and (c), Q2 and x are determined in the phase space which is allowed kinematically. However, a part of the phase space is not allowed by the conditions, e.g. the angular ordering. For the branching at hand, we check the condition and if it is broken, then we cancel the branching, we substitute the maximum virtuality of the parton with its temporal value, Q2, and we return to step (b). (e) Once the branching of the parton is accepted, the remaining job is to select channels of branching and, if the selected channel is three-body, we also determine the light-cone fraction of the second daughter parton. The step is carried out as follows: When the branching parton is a quark, there are three channels, q + g, q + g + g and q + q’ + 4’. Each channel is selected in accordance with the ratio ~(xi):fJ~~(xi,

x 2) dx2:~J~)q~(xi,x2) dx2.

(2.4)

If a three-body channel, e.g. qgg-channel, is selected, then we determine x2 using V~(x1,x2) for fixed x1. When the branching parton is a gluon and it branches into g + anything in step (b), there are three channels, g + g, g + g + g, and g + q + 4. The procedure is the same as that for a quark. For g q4, we need no selection of channels. Now the branching has finished: The mother parton acquires its virtuality. Daughter partons have definite light-cone fractions and they enter step (a) iteratively. The kinematics for two- and three-body branchings is given in sections 3 and 4. Additional constraints in step (d), the angular ordering and the double-cascade scheme, are described in sections 8 and 9. 2, PT). After evolution, we determine the partons’ triplet Q Since thethe Lorentz invariance is respected at everyfour-momenta step in NLLjet,from no their rescaling norvariables, adjusting(x, is required. The determination of p~and p (and P2i. for q4g) in the center-of-mass frame (CMF) for partons at the primary vertex is based on equations is sections 6 and 7. From these values, p p and p ~ of all remaining partons in CMF are determined by equations in sections 3 and 4. The energy and p~are given by eq. (A.2) and the rotational angle ~ is uniformly selected at random. As the last procedure, we rotate the system of partons according to the coupling of the primary quark to the electro-weak current since we develop the parton shower taking the z-axis as the direction of IMF. —*

~,

3. Two-body branching We consider the branching,

p(x= 1, Q2, 0) —~p

1(x,Q~PT)

+p2(1 —x,

Q~,

PT).

-

72

K. Kato, I. Munehisa / NLLjeI: a MC code fore ~e — QCD jets

Here we normalize the light-cone fractions by that of the mother parton p. In general, the transverse momentum of p is not zero. We can write the transverse momenta of Pi and P2 as PT,

XPT + kT,

=

PT,

=

(1

—

x)PT

—

kT.

(3.1)

This equation is used only when 4-momenta of partons are reconstructed, but it is not needed for the argument in the following. With this understanding, the relative transverse momentum, kT, is written as PT in the rest of this section. The squared transverse momentum is given by 2

2

PT~(1~)~Q

I/

2

ç~ ~1

__~__1

The phase spaces for Pt and x~>x>x,

P2

,-~2 ~2

~

(3.2)

are restricted by this equation. The allowed region of x is as follows:

x÷=~(1+e

2.

(3.3)

1—e2±~), D=1—2(e1+e2)+(e1—e2) Here D> 0 is the condition for the allowed region: q~q+g, g—~g+g,

1~2(q+(o)+(cq~to)2>0,

g—q+4,

‘~4~q>°’

2P1P2=

Q~+

1—4~>0,

(3.4)

and

1XQ~+(~T),

(P1+P2)~+1x+X(1_X)~ ~2

~1

Y12~+

or

K

1_x+x(1_x)

Allowed maximum virtuality is

Q~max~(Q2_

~

(3.6)

Q~max(1)(Q2~)

When the virtuality of the sister parton is not yet determined, it is taken to be allowed maximum virtuality.

Q~in order

to calculate the

4. Three-body branching We consider the branching, 2, 0) —~p p(x 1, Q 5(x1,Q~PTI) +p2(~2, Q~PT2) =

+p3(x3, Q~,PT~),

where x1 + x2 + x3 1 and PT, + PT, + PT, 0. Here we normalize the light-cone fractions by that of the mother parton p. In general, the transverse momentum of p is not zero. In this case we can write the transverse momenta of p~(i 1, 2, 3) as =

=

=

PTiPT~T’

(4.1)

/

K Kato, T. Munehisa

NLLjeZ: a MC code fore ~e — QCDjets

73

where k-1.1 + kT2 + kT3 0. This equation is used only when 4-momenta of partons are reconstructed but it is not needed for the argument in the following. With this understanding, the relative transverse momentum, kT, is written as p Ti in the rest of this section. We have =

Q~+PT2 (4.2)

lxi

Allowed maximum virtuality is 2 2 2 ()2 ~J ç~ ~k

Q

Qimaxxi

I

k).

(t~j,

—-~——~-—

k

(4.3)

When the virtuality of the sister parton is not yet determined, it is taken to be Q~to calculate the allowed maximum virtuality. On the other hand, when x 1 is given, the allowed range of x2 is given by x2 =y(1 —x1), —

x3

~Iii5

1

1— x1 —x2,

=

D

—

1

y~>y>y_, (4.4)

4XiQ~min

2—Qfl~ In eq.Y±T (4.4), we± use the fact that the(1—x1)(x~Q minimum virtuality of branching is possible when the following condition holds: —

x

2>Q~ and 1Q q—*q+g+g,

—

P2

is equal to that of

p~.The three-body

D>0,

(4.5)

Q~min’Q~,

Q~minQ~,

q_*q+4+q,

Q~minQ~,

Q~minQ~,

g—÷g+g+g,

Qt~min’Q~,

Q~minQ~,

g—g+4+q,

Qi2min=Q~,

Q~min’Q

(46)

1~

In the case of two-body branching, the transverse momenta are determined by virtuality and light-cone fraction (eq. (3.2)), thereby a trivial freedom for the rotation around the z-axis is left. For three-body branching, we only have eq. (4.2) and PT, + PT, + PT, 0. We determine PT, at random by the uniform 2pT,dp T distribution in the phase space, ~(PT1 + PT2 + PT, )d 2d PT3~The phase space is proportional to =

d~~1 d~~1-2 d~ d~2.

(4.7)

We have the following relations: PT,

=

(x2 +

x3)~T~,

Q~

— =

—

x1

=

PT,

x2+x3 x1

=

~X2PTI

p~-+

_______

x2x

We generate p~,~ and 42 at random for 0 and then we determine PT, with eq. (4.8).

PT,

+PT2,

=

(4.8)

PT2,

~X3~TI

p~-.

(4.9)

2


2 and determine p?~.

2 using eq. (4.9),

5. Primary vertex The term “primary vertex” stands for the process that initiates the parton shower in the e e -annihilation. In NLL order, it is either q4 or q4g. Each parton at the primary vertex has the potential to cascade. +

74

/

K. Kato, I Munehisa

NLLjeI: a MC code fore ~e

—

QCDjets

The multiplicity of partons in a generated event has nothing to do with the primary vertex. In general, the event which starts from q4 has large thrust and that from q4g has small thrust. The mixing of q4 and q4g can be determined in the next-to-leading order. There is a parameter 6 which separates these two types of primary vertex. Though intermediate quantities depend on 6, observables should not have strong dependence on 6 for its moderate variation. This is because we consider all terms in NLL. We write the differential cross section for the ee—8 q4g as da(z1, z2)/dz1dz2, where z~and are 22

the energy fractions of the quark and the antiquark, respectively. We define the phase space S( 6) where the invariant mass of q + g and of ~+ g is greater than 6s. S~=(1—6>z1>0

and

(5.1)

1—6>z7>0).

The above equation assumes that partons are massless. To avoid double counting of the phase space, the 6s. The q~and virtualities of quark and antiquark in the q4 primary vertex are restricted to be less than q4g primary vertex is generated in proportion to the corresponding cross sections. a(q4)=a(total)—a(q4g),

~(q

dz1 dz7.

4g)=fd21~22)

(5.2)

The part of the q4g cross section outside S(6), in which the mass-singularity exists, is absorbed into the next-to-leading order of the Altarelli—Parisi function, P~. The ratio a(q4): a(q4g) and the Pqq (or Sudakov factor which is the integral of P44) depend on 6, but the sum of both contribution not 1~does 2y~//3~ depend on 6 at of least the inclusive limit. This isInthe well-known fact thatwindow the combination C~large. If we is independent theinrenormalization scheme. practice, the allowed for 6 is not take 6 to be very small, a(q4g) dominates over a(q4), which breaks the perturbation. If we take 6 to be very large, most of the generated events start from q4 and the generated distribution is not good in the low —

thrust region.

6. qq primary vertex We consider e~e—~ q(p 1) + P =

-~

=

4i +

1

—

~2

+

4(P2).

\/1

—

The infinite momentum is determined as follows:

2(~+

~2)

+ (e1

—

(6.1)

e2)

It should be noted that I > > 0 and that s~ I when Q~ 0. The inside of the square root in eq. (6.1) is positive in the region bounded by a parabola that passes through (e’, e2) (0, 1), (~,~),(1, 0). The components of four-momenta Pi and P2 are =

=

=

=

(1

—

~)W.

p~=(e1/s~)W,

PT°’

p~ (I

PT,

=

—

c1/~)W,

=

(6.2)

As discussed in the previous section, there is the restriction on virtuality as

6s> If 6

<

Q~,Q~.

~,the inside

of the square root of eq. (6.1) is positive.

(6.3)

K Kato, I. Munehisa

/

NLLjet.’ a MC code fore ~e — QCDjets

75

7. q4g primary vertex We consider e~e—*q(p1) + 4(P2) + g(p). The transverse momenta of Pi and p3 are The infinite momentum is determined as follows: =

4i +~

=

(Pt YY13

2(Y + ~2) —

~2 2

+

+~3)

~1

—

~PT•

~2)

(7.1)

K

(3

and

s

±

=x~W,

1+K

Pi

(1 —~)W,

> 0 and that sj

=

1 when

Q~ =

0. The four-momenta of Pi’

P2

(

—

=

(~

—

It should be noted that 1 > y> 0 and 1 > and p~are

Pt

+

PT

W,

=

p~ (i =

PT.

(7.2)

PT2°,

—

(3

PT’

+ K

p~=(1—x1)~W, ~P3 (l—x1)~”

PT•

PT3

When three partons are massless, we have Q~=Q~=Q~=0, p~=xW,

p~=(1—x)W,

i~=1,

pj =t(1 —x)W,

PT,

z1 =x+t(1 —x),

p~=(1—t)W,

PT2°,

z2=(1—t),

p~=txW,

PT3= ~PT’

z=1—x+tx,

PT’

(7.3)

where x=x1

and

(7.4)

t=(1).

Other relations for the massless case are z1 + z2

+

z3

=

2,

Y13

=

1

=

t,

Y23

=

1

—

=

(1— x)(1

—

t).

(7.5)

There is a delicate problem in connection with the q4g primary vertex. To generate this vertex we use the cross section dG(z1, z2) which is calculated for the on-mass-shell q, 4, and g. Instead of z5 and z2, we can use x and t and determine them from do( x, t). However, these three quantahave acquire virtuality 2~PT), the three partons (x, 0, PT), (1,in0,the 0) course their cascade.The In antiquark the triplet cascades variables,backward (x, Q and (1 of x, 0, —PT). to the quark and the role of p~ and p is reversed for the cascade on the antiquark side. Later, the virtualities will get non-zero values. Consequently, we have an ambiguity when we calculate their four-momenta though it is of the nature of the next-next-to-leading order. We have two ways of doing this. —

(a) The invariant mass kept constant. From x and PT we calculate Y13 and Y23 by eq. (7.5). From s, Q~,Q~,Q~,Y13 and Y23’ we calculate 4-momenta assuming that y is defined by eq. (A.9). The transverse momentum of P2 ‘5 0 and p~is

76

K. Kato, I Munehisa

/

NLLjet.’ a MC code fore e

-

QCD jets

replaced by Pt:

(PIP2

=

—

E1E2)/I P21’

P3:

(P2P3

=

—

E2E3)/I

P2’

PT—

I

Pi:

=

P~I (7.6)

(b) The angle kept constant. From x and PT we calculate the energy fraction z~using eq. (7.3). We assume that the momentum fraction is equal to the energy fraction. We use the relation for the massless case to determine angles between partons, i.e. cos 0, I 2(z1 + z, 1)/(z7z1). Next, solving the equation W= ~(Q~+ on a, we get the absolute value of space momentum of p, p I ~/~z1W/2. From the above cos 9~,and I I we can calculate 4-momenta assuming that the transverse momentum of P2 of 0. =

—

—

=

We use method (a) in NLLjet as a default. However, method (b) and a “mixed” method of (a) and (b) are available as options. Finally it should be noted that at every step a check must be made to ascertain that the inside of the square root in eq. (7.1) is positive.

8. Angular ordering The angular ordering in our scheme is fully described in ref. [10]. It is the higher-order effect and we can show it explicitly through the study in NLL. In some codes for the parton shower, the angular ordering is introduced for every branching. As shown in ref. [10] this is not correct. Our angular ordering originates from a careful analysis of the three-body vertices, (a) g(p) g(p1) + g(p2) + g(p3), and (b) q(p) q(p1) + g( P2) + g( ps). The light-cone fraction is normalized so that x 1 for p. For quark branching, there is no restriction due to angular ordering. We consider a gluon of momentum P23 which branches into P2 +p3. When the mother of the gluon P23 is a gluon (p), the vertex (a) is applied and when it is a quark (p), then the vertex (b) is applied. As is described in section 2, the condition is applied after we determine the virtuality of P23’ Q~3.and the light-cone fraction of P2 and p.~,x7 and x1, in the conventional way. --~

—~

=

(a) Mother parton is a gluon: >

Q~< Q~X1X2

x2, x3,

x1x2x.1 X2X3

~

X3X1

(1

+

~R log R),

R

x2x3 =

-

2

(x~+ x~)

(8.1)

It should be noted that when the first inequality does not hold, there is no restriction. (b) Mother parton is a quark:

Q~3<

x2x3 X1X2

+

X2X3

+

X3X1

(1

+

7R log R),

R

=

x-,x1 (x + x

)2

.

(8.2)

9. Double-cascade scheme and the separation of jets As described in ref. [11], we use the “double-cascade scheme” to improve the symmetry between quark jets and antiquark jets. The scheme comes from the observation that eqs. (7.3)—(7.5) are symmetric under

K. Kato, I Munehisa

/

NLLjeI: a MC code fore ~e

—

QCDjets

77

the interchange of t and 1 x when I is small. Phase space parameterized by x and t is divided into 1 x> t and 1 x < t. The first part corresponds to jets on the quark side and the second one corresponds to jets on the antiquark side. Using the symmetry, the role of 1 x and t is interchanged in antiquark jets. Then —

—

—

—

(9.1)

t<1—x

is required at every branching. This condition works for separating partons in a quark jet and those in an antiquark jet in an event that originates from primary q4. For primary q4g-jets there is no explicit condition; constraint on the phase space in q4g-jets is effective only beyond the next-to-leading order. However, the introduction of a condition to separate each jet seems to be natural. Hence, we assume that eq. (9.1) or an other equivalent expression should hold at every branching. The practical implementation of eq. (9.1) is done as follows. We consider partons of momentum k, p~and so on. The parton to be tested by eq. (9.1) is the mother parton of parton k, i.e. we test the branching which produces parton k and its sister partons. We assume that p, is inside of the jet which contains k and that p~is outside of it. There are many equations that are equivalent to eq. (9.1) in the massless limit; p,,

(Al)

kp~
(A2) (Bi)

(k+p1)2<(k+p~)2, kp 1/E~
(B2)

1—cos9<1—cos91,

(92)

where (Al) and (A2) are the condition for separation by invariant mass, and (Bl) and (B2) are that by angle. We call these equations the condition of the angular separation (of partons). When k, p, and p~are massless, (Al) and (A2) are the same and (Bl) and (B2) are the same. In the massless case, we call them simply (A) and (B). 9.1. q~primary jets, direct branching The “direct branching” is defined as the branching of a parton that+belongs to vertex. The We 2, PT)) g(p (1 thex,primary Q~’—PT)). consider e~e—~ 4(P2) and way. q(p1) q(k parton (x, K is p branching of P2 is q(p~)+ treated in the same The inside 3 or Pi and the outside parton is P2~By using of the equations in section 6 we have =

—*

2kp3 s

l—x x

+

l—x

sc

-

(3+

~-)(1—

~.i~2xi~(l_ x

x l—x

—

=

x(l—x)

i,) ~

(9.3)

K

x(1—x)

The last one comes from eq. (3.5). The interchange of quark (k) and gluon (ps) gives us another equation obtained by replacing x with 1 x.2 Q~ 0, we have If we take them massless, i.e. K (A),(B)~ x(l~x)
=

=

=

78

K. Kato, I Munehisa

/

NLLjeI: a MC code for e + e — QCDjets

For x 1, the second relation reduces to eq. (9.1), which is identical to the condition to separate the gluon from the antiquark in this limit. Since (1 x)/(2 x) < 1 x, eq. (9.4) is a stronger condition than eq. (9.1). —

—

—

9.2. qijt primary jets, general We consider a parton k (xxm, K2, kT) that belongs to the cascade of the quark (Pi). Here Xm is the light-cone fraction of the k ‘s mother-p arton and in the rest of this section x is the relative light-cone fraction. The sister parton of k, k’ ((1 x)xm, K’2, k~-),is also considered implicitly as we exchange x and 1 x. We consider the separation with respect to Pi (in) and P2 (out). =

=

—

—

2kp 1

=xx m + ~

Here ~cis for

kT.

~k

+

2kp2

K

XXm

,

ii

+

K

XXm71

(9.5)

.

In order to compare the results in subsection 9.1, we take

lxxx~i(12dl)<

~m<

—) + (I ~) 1

/

=xxmii~l

—

S

k

=

~

and

ij

I for case (A):

=

l—~’

(9.6)

where m x(l x)K is the normalized virtuality of the mother parton. As x 0, eq. (9.6) reduces to eq. (9.1). We get a similar relation for a parton k that belongs to the antiquark cascade after the exchange of 1 x and t. =

—

—

9.3. q~gprimary jets, general 2, kT) (and its sister parton). We we consider q(p1)the + 4(P2) + g(p3) and a7,parton k Again take e~e-~ k 0. Using equations in section we have =

2kp 2

5~Xq

(~+ K

—

xq

Xm

2 I x —-=(l—~)~ Xm —

s

m

(xxm, K

2kp

2kTpT 1

=

‘

/

y\ (9.7)

where Xq is the light-cone fraction of the quark at the primary vertex and K is the normalized k~.We neglect the last term of the first equation of eq. (9.7). The parton k belongs to the quark or antiquark cascade. The condition for the separation of k from p2 or Pi is given by k in quark,

X (

<

1

—

2 1 x x mxq 1 —

2e~ —

x1

(2_Xq)K —

xq(1 _xq)

__________

(9.8) x

kinantiquark,

Cm>

lxxrn

K(2—xq) 1+ flxq(l_xq))

Using these relations we define s’, the effective s, as is shown below. 9.4. Summary We summarize the results of this section. The condition of the double-cascade scheme is 2/s)<1-x (Q

(9.9)

K Kato, I Munehisa

/

NLLjet: a MC code fore ~e — QCDjets

79

for the branching q(x 1, Q2) q(x) + anything. If the event starts from the q4g primary vertex, we replace s in eq. (9.9) by s’ where s’ is given by —~

=

from quark

s

=

s

from antiquark

s’

=

s 1

Xq

+

PT

(9.10)

2 —

Pr

xq(l _xq)s

The angular separation eq. (9.2) is a stronger condition than the condition imposed by the double-cascade scheme. We do not apply eq. (9.9) for the branching of the gluon since the angular-ordering condition in the section 8 is imposed for gluons. For the branching of partons that belong to the primary vertex, we further impose the angular separation to establish the separation of jets. This is because the double-cascade scheme is not a strong enough condition and it does not work very well for massive quarks or highly off-mass-shall case. For the q4-case we require that

IkHIp

2I

(9.11)

in the notation of subsection 9.1 and the q4g-case we require (Bl) in eq. (9.2).

10. Determination of momenta

2~PT). The As described in section 2, weand generate four-momenta from triplet variables (x,formulation Q direction of p defines the IMF is related to the selection of their the gauge vector. In our the vector of the light-cone gauge is parallel to the momentum of the primary antiquark and the antiquark does not radiate a gluon. The introduction of the double-cascade scheme and the existence of hard three-jets change the situation. Relations in sections 3 and 4 hold for a specific direction of p Using these relations, we first determine p and PT’ and then determine p E and p~are determined by p and p. The four-momenta of partons are determined as if p~points the positive z-axis. After that, we rotate partons according to the direction of the primary parton. The polar angle of a primary parton is defined as follows: +

~.

+

.

+

(a) q4 9=0;

4:

9=rr;

q:

9=0;

4:

9=’rr;

g:

0

q:

(10.1)

(b) q4g

for

(Pi +p3)

<(P2

+p3),

IT

for (Pt +p3)> (P2 +p3),

(10.2)

where Pt’ P2 and p~are the momentum of the quark, of the antiquark, and of the gluon, respectively. For the q4-case, (a) is the only possible choice. However, for the q4g-case, we can define a different scheme (b’). (b’) q4g polar angle of P~. (10.3) In NLLjet the default method is (b’) but (b) is also available as an option. It is noted that the difference between (b) and (b’) is beyond the NLL order. =

80

K. Kato, I Munehisa

/

NLLjet: a MC code fore ~e — QCDjets

11. Mass of the quark Proper treatment of the quark mass, m, is beyond the validity of the present formulation since it is constructed for massless quarks. However, we must, somehow, take it into account. For instance, the probability of the branching of a gluon into a quark pair must be flavor-dependent since the string picture suggests that the probability of the pair creation of quarks from the vacuum is suppressed as exp( cm2). Moreover, if we neglect the quark mass, we may have a trouble with their hadronization as a mismatch between the mass of the produced hadron and that calculated by partons’ momenta. In order to cure these defects, we have introduced m into NLLjet. We list below the places where m appears. Some of them are very natural conditions, some are essential in order not to break kinematics, but some are only empirical relations and need further elaboration. —

(a) Primary vertex: We reject a configuration that is not allowed kinematically when we turn on the quark mass. (b) Modification of

Q~: (11.1)

Q~=Q~+m2.

Since the minimum of Q2 is changed, it affects the evolution of the parton shower in many steps. For instance, the probability of the branching g q4 depends on the flavor of the quark since the available phase space is dependent on m. -~

(c) Modification of 6: 6=maxof(6,2~~/W).

(11.2)

When the radiative correction is turned on, this relation must be tested in each event since W is not constant. (d) Double cascade: We modify condition (9.9) as (Q2/s)

-

(m/W) <1- x,

where W = ~ for q4 and

W=

(11.3)

v~i~ for q4g.

12. eke-physics The coupling of a quark to the Z°-gaugeboson is defined as j(Z)

=

2~/x~(l— x~)4[I~y,,(l

—

y 5)

—

2eqxwy~~q,

(12.1)

20~)is the electro-weak mixing angle, 13 (= ±~) is the weak-isospin of the quark, and eq where x~(= of sinthe quark in units of positrons. The mass of Z°,the width of Z° and x~are independent is the charge parameters in NLLjet and they can be specified by the input parameters. The angular distribution of the primary partons is determined using the lowest-order martix elements, i.e. e~e—~ y, Z0 q4 for the primary q4-vertex, and e~e—’y, Z° q4g for the primary q4g-vertex. The production rate of each flavor is calculated by the former. —*

—*

/ NLLjet:

K. Kato, T. Munehisa

a MC code fore ~e — QCDjets

81

The radiative correction to the e e -state is important. NLLjet includes a function to generate the initial-state radiation based on a formula similar to that of Berends, Kleiss and Jadach [19].The final-state radiation is neglected here. The modules for the emission of a hard photon are designed on the assumption +

that the energy dependence of the hadronic cross section is moderate. This assumption is allowed in the PEP/PETRA/TRISTAN region. In the SLC/LEP region, however, the generation of the hard photon is less efficient so that it is recommended to switch off the flag for the radiative correction. Input parameters are open to users as explained in appendix D. We have checked that the physical results do not strongly depend on 6 nor Q0. However, the large deviation of these values from their default values should be avoided. For example, we face non-negligible model dependence on hadronization when we take Q0 very large; we may break the perturbation when we take Q0 or 6 very small; we get wrong distribution in large PT when we take 6 very large. For a standard, we suggest 0.08 <6 <0.20 and 0.5 < Q0 < 4.0. In order to study experimental data in e e -reactions and to determine the QCD A, a hadronization model must be linked with NLLjet. Analysis of MARKII data at W = 29 GeV and other sets of data was done by Kamae and Shirahashi using the LUND string fragmentation [16]. They get A ~ = 235 ±52 MeV. Analysis of VENUS data at TRISTAN was also done using EPOCS hadronization [17] and a value of A~= 254i~±56 (sys.) MeV was obtained [18]. +

Acknowledgements

The authors are grateful to Y. Shimizu, T. Kamae, A. Shirahashi and T. Watanabe for useful discussions and for their help with the development of the program. One of the authors (T.M.) would like to thank T. Sjöstrand for his advice and the CERN Theory Division for its hospitality. This work is partly supported by the Special Research Fund of the Kogakuin University, Japan.

Appendix A. Basic notations We use two coordinate frames. The first is the infinite momentum frame (IMF) and the second is the center-of-mass frame (CMF). The evolution of the parton shower is executed in IMF. The direction of the boost points to the positive z-axis and the transverse plane is the x—y-plane. 2= = squared c.m. energy (A.1) s= q p~=E+p~, p=E—p~, E=~(p~+p), p,=~(p~—p). (A.2) Q2 =

virtuality,

Q~ =

minimum virtuality or minimum virtuality of the gluon,

Q~ =

minimum virtuality of the quark. Q2+p2.

-

kT,

x

PT

=

transverse momentum,

P~PT’

~2

=

4-

=

p

=

p •p’ = ~(p~p’+pp’~)

light-cone fraction in IMF,

(A.3)

.

~PTP’T•

P = infinite momentum.

(A.4) (A.5) (A.6)

K. Kato, I Munehisa

82

/

NLLjet: a MC code fore ~e

QCDjets

We define some dimensionless quantities: (A.7) that are used as e0, or

Cq~ k

for

Q~,

Q~,

K~,respectively, and as e~for a momentum p,. (A.8)

~,

that may have subscripts similar to those of c. Also ~~iJ=PT

PT/S,

and similar notations are used.

y=(Pi±Pi)

=

~

=

energy fraction in CMF.

(A.9)

(A.I0)

Finally, we define color factors for Ncoior = 3. CF=~, 3, CA= TR=~,

TaTa=CF•I,

fahcfahdCA6cd,

(A.I1)

Tr(7~Tb)=TR6Oh.

Appendix B. List of subroutines In this section, we present the modules in NLLjet with a short explanation of each module. The first line shows the SUBROUTINE statement or the FUNCTION statement that just appears as the first (non-comment) statement of the module in the NLLjet program. For an argument, we show its name in capitals and sometimes add a tag (in) and/or (out) to specify that the argument is referred and/or defined in the module. Modules and variables used in the program are also shown in capitals. Before the generation of events, two routines, QCINPT and QCPREP, are to be called for the input of parameters and the initialization. After that, QCGEN is called repeatedly to generate events. Samples of the main program are found in appendix E. SUBROUTINE QCINPT(WBEAM, NEVEN, NRADC, NPRNT) Reads the control parameters for NLLjet from the file assigned to unit FTOSFOO1. WBEAM(out) is the center-of-mass energy, in GeV, NEVEN(out) is the number of events to be generated, NRADC(out) is a flag for the radiative correction (NRADC = IRC in input parameters) and NPRNT( out) is the interval to print the list of generated events. BLOCK DATA QCDATA Sets the default values of parameters that are to be read in QCINPT. SUBROUTINE QCPREP Controls the prologue procedure for NLLjet. It calculates the variables in common areas and calls the modules for calculating the table of the Sudakov factors and the cross sections for the primary vertex. In order to save computer time, the table and the cross sections may be read from a file if they have already been calculated and saved in a file.

K. Kato, T Munehisa

/

NLLJet: a MC code fore~e - QCDjets

83

SUBROUTINE QCOUT(NEVE, NCHK, WCMS) Prints the four-vectors of partons (PARTN) and the status of partons (LISTP) in the generated event. If NCHK(in) = 1, the momentum conservation is examined assuming that the event is generated in CMF with the total energy of WCMS( in). NEVE( in) is the sequential number of the generated event. SUBROUTINE QCROT The parton shower evolves taking the z-axis as the direction of the IMF. In this routine, four-momenta of generated partons are rotated using the distribution calculated by the coupling of the quark to the electro-weak current. SUBROUTINE TBFMD(RS, QM, GM, NWEAKX) Calculates the parameters for the angular distribution of two-jets. These parameters are used in QCROT. REAL FUNCTION FJET2(COST, NF, ITOT) This function gives the differential or the total cross section of two-jets. NF specifies the flavor. When ITOT =0, FJET2 = (l/~)da/dQ for cos 9 = COST and when ITOT = 1, FJET2 = (1/a~)a. REAL FUNCTION FJET3(X, NF, NOK) This function gives the differential cross section of three-jets. NF specifies the flavor. FJET3 z1)(1 z2)do/dz1dz2dQ1d42 in unit of

=

(1

—

—

SUBROUTINE TBPRV This subroutine calculates the fraction of primary two-jets, S2MOD in /PRVTX/. SUBROUTINE JR31ST(DELT, GO) Calculates the fraction of primary three-jets of order a8 by integrating the differential cross section in leading order. SUBROUTINE JR32ND(DELTA, NL, Gi, GR1) Calculates the fraction of primary three-jets of order a~and the average of ~ numbers that have been calculated numerically for some 6 = DELTA( in).

.Y2

by interpolating the

SUBROUTINE INITA This subroutine is responsible for the2).It following initialization. It prepares coefficients forbythethe calculation prepares the array for numerical integration Simpson of the strong constant, a~(Q method that iscoupling used in PNExx to calculate the Sudakov-factor tables. It also checks the mass of the quarks in order to save CPU time. SUBROUTINE TBOUT This is a service routine to print the Sudakov-factor tables. SUBROUTINE TBPNE The table of non-branching probability, i.e. the Sudakov factor, is calculated in this subroutine using subroutines PNExx that are explained in the following. The tables are stored in /TPNExx/. SUBROUTINE PNE1A(QQO, QQ2, NF, PNE) Calculates the leading part of the Sudakov factor for q = QQ2) is done analytically.

—~

q

+

X. Q2-integration from

Q~( =

QQO) to

84

K. Kato, I Munehisa

/

NLLjet: a MC code fore ~e - QCDjets

SUBROUTINE PNE1B(QQ1, QQ2, NF, PNE) Calculates the sub-leading part of the Sudakov factor for q Q~(= QQ1) to Q~(= QQ2). SUBROUTINE PNE2A(QQO, QQ2, PNE) Calculates the leading part of the Sudakov factor for g Q2( = QQ2) is done analytically.

-~

g

—~

SUBROUTINE PNE2B(QQ1, QQ2, PNE) Calculates the sub-leading part of the Sudakov factor for g Q~(= QQ1) to Q~(= QQ2).

+

q

+

X by numerical integration from

X. Q2-integration from

-~

g

+

Q~( =

QQO) to

X by numerical integration from

SUBROUTINE PNE3C(QQ1, QQ2, NFLAV, PNE) Calculates the Sudakov factor for g q + 4 by numerical integration from Q~(= QQ1) to Q~( = QQ2). —~

SUBROUTINE PNELL

The above subroutines PNExx use the fact that the argument of a 8 is p~= x(l in the leading order for a 2). 8(Q

the Sudakov factor

2. This calculates —

x)Q

SUBROUTINE QCGEN(IFLAV) This subroutine controls the generation of one event. The meaning of IFLAV(in) is explained in GENPRV. First, it resets the common area /QCEVNT/ that records the generated event and it calls GENPRV to create primary partons. It iteratively calls GENQBR and GENGBR that process the branching of the quark and of the gluon, respectively. After the evolution ends, the four-momenta of partons are determined: It calls CMSPRI in order to register ancestors and calculates the four-momenta of partons at the primary vertex by CMSPRV. After the directions of primary partons are recorded in /CMSANG/, the partons are forced to forward the positive z-axis. CMS2BD and CMS3BD are called iteratively to determine four-momenta of partons at each two- and three-body branching. The direction of partons is recovered by CMSROT. Finally, the event is rotated by QCROT according to the coupling of the quark to the electro-weak current. SUBROUTINE GENPRV(IFLAV) Determines the type of primary vertex and generates it for either the q4-state or the q4g-state. If IFLAV(in) is zero, the quark flavor is selected according to the R-ratio. Otherwise, events with flavor = IFLAV are always generated. SUBROUTINE GENQBR(NQUA) NQUA( in) is the address of the quark to be processed in this subroutine. The quark is tested for branching or no-branching. If it branches, then its virtuality, the mode of branching, and light-cone fractions (x) of daughter partons are determined here. SUBROUTINE GENGBR(NGLU) NGLU( in) is the address of the gluon to be processed in this subroutine. The gluon is tested for branching or no-branching. If it branches, then its virtuality, the mode of branching, and light-cone fractions (x) of daughter partons are determined here. SUBROUTINE PRVQQG(X, PT, NF, MCND) Determines the x(= X(out)) and I PT I(= PT(out)) for primary three-jets of flavor = MCND = 1, this means that the trial in this routine has failed.

NF(in).

When

K Kato, I Munehisa

/ NLLjeI:

a MC code fore ~e

—

QC’Djets

85

SUBROUTINE CMSPRV(NERROR) Determines the four-vectors of the primary partons from given squared virtual mass and fractions x of light-cone variables. It returns with NERROR = 1 if the kinematics for three-partons is broken. After recording the direction of partons (section 10) in /CMSANG/, all primary partons are rotated to evolve along the positive z-axis. See sections 6 and 7. SUBROUTINE CMS2BD(NADR) Determines the four-vectors of the daughter partons in the two-body branching from the four-vectors of their mother parton and virtualities and x of daughter partons. NADR( in) is the address of the mother parton. The azimuthal angle is chosen at random. See section 3.

SUBROUTINE CMS3BD(NADR) Determines the four-vectors of the daughter partons in the three-body branching from the four-vectors of their mother parton, and virtualities and x of daughter partons. NADR( in) is the address of the mother parton. Variables that cannot be fixed, namely, azimuthal angles, are chosen at random. See section 4. SUBROUTINE CMSPRI By using LISTP(5, ), this subroutine registers the ancestor of each parton. The ancestor is the parton at the primary vertex which produced the parton in question. Thus LISTP(5, ) = 2 or 3 (q or 4~for the event that starts from the two-jet primary vertex and LISTP(5, ) = 2, 3 or 4 (q, 4 or g) for that which starts from the three-jet primary vertex. SUBROUTINE CMSROT According to /CMSANG/, the direction of the partons is restored in this subroutine. SUBROUTINE LIMAO(X3, RR, Z, MOTHER, NOK) This subroutine checks whether the branching is valid with respect to the condition of angular ordering. NOK(out) = 1(0) means that it is accepted (not accepted). See section 8.

SUBROUTINE LIMDC(NQUA, NID, K2, Z, NBR) This subroutine checks whether the branching is valid with respect to the condition for the doublecascade scheme. NBR( out) = 1(0) means that it is accepted (not accepted). See section 9. SUBROUTINE LIM3(KCASE, NADD, NBR, Q2, X) This subroutine checks whether the generated primary three-jet system is consistent with the kinematical condition. NBR(out) = 1(0) means that it is accepted (not accepted). When KCASE(in) = 0, the test is done for the three partons that belong to the primary vertex. When KCASE = 1, the test is done for the branching of one of the three partons in the primary vertex. See section 9.

SUBROUTINE LIM2(NADD, Q2, Z, NBR) This subroutine checks whether the generated primary two-jet system is consistent with the kinematical condition. The test is done for the branching of a parton at the two-jet primary vertex. See section 9. SUBROUTINE LIMBR(NADD, Q2MAX) Determines the maximum virtuality Q2MAX( out) for the branching of a parton at address NADD( in). See eqs. (3.6), (4.3) and (6.3).

K. Kato, I. Munehisa

86

/

NLLjet: a MC code fore ~e -. QCDjets

SUBROUTINE CHKKN1(S, Y12, Y23, Q12, Q22, Q32, NCND) Tests whether the three-parton state is kinematically allowed. NCND(out) is satisfied. See item (a) in section 7.

=

0 tells that the kinematics

REAL FUNCTION EFF3JT(YA, YB. YCUT, ALPHA, NL, NSL) This function gives the differential cross section of three-jets up to the order ~ REAL FUNCTION EF3JTF(YA, YB, YCUTO, ALPHA) This function is called from EFF3JT and gives the non-singular contribution of the order a~differential cross section of three-jets. FUNCTION EF3JCT, EF3JAO, EF3JH, EF3JK, EF3JQO, EF3JFC, EF3JFN These 7 functions are parts of EF3JTF. REAL FUNCTION APQX(Z, K2) The Altarelli—Parisi function for the inclusive decay of a quark of light-cone fraction x( = Z( in)) and virtuality Q2(= K2(in)). REAL FUNCTION APQXI(Z) This function gives the value of the second-order Altarelli—Parisi function of the quark. REAL FUNCTION APQX1A(Z) This function gives the leading term of APQX1 around x = 1. REAL FUNCTION APGX(Z, K2) The Altarelli—Parisi function for the inclusive decay of a gluon of light-cone fraction x( = Z(in)) and virtuality Q2(= K2(in)). REAL FUNCTION APGQQ(Z, K2) The Altarelli-Parisi function for g Q2(= K2(/n)).

-~

q4 for a gluon of light-cone fraction x(= Z(in)) and virtuality

SUBROUTINE HAPQX(K2, ZMIN, ZMAX, ALPH, Z) Determined the light-cone fraction x ( = Z( out)) for ZMIN( in)
1

—

EPSL at the branching of a

for EPSL( in)
FUNCTION V3QGG(Z, Y), V3QQQ(Z, Y, LF), V3GGG(Z, Y), V3GQQ(Z, Y) These functions give the probability distribution of the three-body branchings: q

—*

q(x

1) + g(x7) + g(1 q~q(x1)+q’(x2)+~’(1 —x1 —x2), g~g(x1)+g(x2)+g(1 —x1 —x2), and g~g(x1)+ x1 x2), respectively. Here, Z(in) = x1 and Y(in) = x2. For V3QQQ, LF(in) = 1(2) means that q ~ q’ (q = q’).

—x1 —x2), q(x2) + 4(1

—

—

K Kato, I Munehisa

/

NLLjet: a MC code fore ~e — QCD jets

87

FUNCTION SV3QGG(Z, EPSL, ALPH), SV3QQQ(Z, EPSL, ALPH, LF), SV3GGG(Z, EPSL, ALPH), SV3GQQ(Z, EPSL, ALPH) These functions give the value of the integral of the three-body branching function V3QGG, V3QQQ, V3GGG and V3GQQ, respectively The integral is done by the variable x2 fixing x1( = Z). For SV3QQQ,

the meaning of LF(in) is the same as that in V3QQQ. The argument ALPH(in) is the value of the strong coupling constant and EPSL(in) is unused at present. SUBROUTINE SV3GW, FiN, FJQ, FJY, FDF, FMF These six functions are the parts of SV3GGG.

SUBROUTINE HV3QGG(Z, YMIN, Y), HV3QQQ(Z, YMIN, Y, LF), HV3GGG(Z, YMIN, Y), HV3GQQ(Z, YMIN, Y) These subroutines determine the light-cone fraction of the second daughter parton, x2( = Y( out )), for the three-body branching q q(x1) + g(x2) + g(l x1 x2), q q(x1) + q’(x2) + 4’(l x1 x2), g g(x1) + g(x2) + g(l x1 x2), and g g(x1) + q(x2) + 4(1 x1 x2), respectively. The light-cone fraction of the first daughter parton, x1(= Z(in)), is fixed here. For HV3QQQ, the meaning of LF(in) is the same as that in V3QQQ. —~

—

—

—

—~

—

—

—‘

—

—

—*

—

REAL FUNCTION ALPS(QQ, MODE) 2). The argument QQ(in) is Q2 in GeV2, function the strong andThis MODE( in) =calculates 0(1) specifies that a coupling constant, a8(Q 8 isand calculated using the leading 2/A2), in the next-to-leading case,(next-to-leading) order formula. In the leading case, a,, = 1/c0 log(Q =c

0 log(3-) +c1 log(~_+~1)

is solved iteratively where c0 and c1 are CO and Cl in /COMAL1/. See appendix C. REAL FUNCTION ALPX(QQIN), ALPY(QQIN, Z), ALPD(QQIN) 2) in many These ALPY functions givea~(x(l strong coupling at ALPD = QQIN(in). ALPX iscontribution used as a8(Q routines. is the x)Q2) forconstant x = Z(in). is the non-leading of a~(Q2). —

SUBROUTINE HBR1DM(PNE, TP, TMIN, TMAX, DELT, NTMX, Ti, NBR) Determines virtuality at a branching according to the table of the probabilities PNE(in). The virtuality is given by A2 exp(Ti(out)).

SUBROUTINE INTRPL(PNE, T, TMIN, DELT, NTMX, P) Gives the interpolated value P(out) for T(in) using the table PNE(in). REAL FUNCTION SPENCE(X), SPO(X), COMPLEX FUNCTION CSP(Z) These functions give Spence function, sp(x) = —J~dtlog(i t)/t, for x SPENCE and CSP is complex-valued for the complex argument Z(in). —

=

X(in). SPO is used by

REAL FUNCTION P2BODY(AMOO, AMO1, AMO2) This function gives the momentum of the two-body decay in the center-of-mass system for a particle with a given mass AMOO(in) to particles with masses AMO1(in) and AMO2(in).

88

K. Kato, I Munehisa

/

NLLjet: a MC code fore ~e

-

QCDjets

SUBROUTINE BOOST$(PPPP, PREP, NNN1, NNN2) This subroutine executes the Lorentz boost of a set of four-vectors, PPPP(4, I)(in and out). For a given four-vector PREP(4)( in), the system is transformed into another system in that the space component of the four-vector PREP(4) is 0. Here, I runs from NNNI(in) to NNN2(in). SUBROUTINE ERRPRI This subroutine controls the printing of error messages. It has a number of entry points ERRPxx. SUBROUTINE ROTAT$(PPPP, PHIl, COST, PHI2, NUMP, NDIM) This subroutine executes the Euler rotation of a set of vectors. SUBROUTINE SWAPD(X, Y), SWAPI(X, Y) These subroutines interchange two arguments X(in and out) and Y(in and out). SWAPD is for real variable and SWAPI is for integer variables. REAL FUNCTION RANDOM(IDUM) Generates a uniform random number between 0 and 1. IDUM(in) is a dummy argument. This module is dependent on the system. SUBROUTINE XXDATE(CDATE) To record the run, it provides some information, e.g. the date of the run, as CDATE( out) that is a string of 8 bytes. This module is dependent on the system.

SUBROUTINE QCRCGN, QCRCBS, RCGEN$, RCGIN$, RCGCN$, RCGRS$, RCGAN$, RCICN$, RCIRS$, RCCOE$, RRATF$ These subroutines are responsible for the radiative correction for the initial state radiation. The method is similar to that in ref. [19]. A detailed explanation is skipped here. The usage is shown in appendix E.

Appendix C. Common areas In this appendix we present all common areas used in NLLjet with a short explanation. The name of common areas and that of variables are shown in capitals. /QCEVNT/ PARTN(10, 1000), LISTP(lO, 1000), NATOP, NSIZE Common area to record particles in the generated event as a list. PARTN: PARTN(J, I) records momentum of the parton at address is as follows: PARTN(1, I) x, the fraction of the light-cone variable. PARTN(2, PARTN(3, PARTN(4, PARTN(5, PARTN(6,

1) I) I) I) I)

Q2,

I

=

I. The meaning of the first index, J,

the virtuality, in GeV2.

I~the magnitude of transverse momentum, in GeV. 4, azimuthal angle of PT• p~,x-component of momentum, in GeV. p 8,, y-component of momentum, in GeV. PARTN(7, I) p~,z-component of momentum, in GeV. PARTN(8, I) E, energy, in GeV. PARTN(9, I) p ‘-component of momentum, in GeV. PARTN(1O, I) p-component of momentum, in GeV. PT

K. Kato, T. Munehisa

/ NLLjet:

a MC code fore ~e — QCDjets

89

LISTP: LISTP(J, I) records the history and the status of the parton at address = I. The meaning of the first index, J, is as follows: LISTP(l, I) Code number of parton. 0 is gluon, and 1, 2, 3, 4, 5 and 6 are the u-, d-, s-, c-, b- and t-quark, respectively. For antiquark, the code number has a negative sign. The primary vertex whose address is 1 has a special code number, 99. When the radiative correction is switched on, the radiated photon may be added at the bottom of list with the code number 100. LISTP(2, I) Relative address of the mother-parton. The address of the mother parton is I-LISTP(2, I). LISTP(3, I) Number of the daughter partons. If it is 0 for a parton, the parton is a final-state parton. LISTP(4, I) Relative address of the eldest daughter parton. The addresses of the daughters are from I + LISTP(4, I) to I + LISTP(4, I) + LISTP(3, I) 1. LISTP(5, I) Status of the parton. This is used in two ways. In the phase of parton evolution, it is used as follows: 0 means x and Q2 are undetermined; 1 means x is fixed but Q2 is undetermined; 2 means x and Q2 are fixed. After the evolution, 4-momenta are determined by x and Q2. It is used to record its ancestor parton. —

NATOP: The top address in event. NSIZE: Size of PARTN and LISTP. At present, it is 1000. /NLLPAR/ ROOTSO, ROOTS, Q2MAX, Q02, FLAMD, FLMD2, DELTA(lO), QMASS(iO), SCHEME, NFLAV Important physical parameters for NLLjet. ROOTSO: Input value of the center-of-mass energy, in GeV. ROOTS: Center-of-mass energy, in GeV, of the generated event. It may differ from ROOTSO due to radiative correction. In physical routines of NLLjet, ROOTS is always referred. Q2MAX: Square of ROOTS, in GeV2. Q02: Minimum virtuality, Q~, in GeV2. FLAMD: QCD A parameter, in GeV. FLMD2: Square of FLAMD, in GeV2. DELTA: Parameter 8 that separates the hard scattering and the shower in the first branching. It may depend on the flavor of quark. QMASS: Mass of quark, in GeV. SCHEME: Scheme parameter of the P-function for quark branching. NFLAV: The number of flavors, Nf. /QCNLL/ NLAP, NLAL, KZAP, KZQQG, NSHUF1, NSHUF2, NLMODE, NAO Control parameters for NLLjet. These are flags to select various modes. NLAP: Selection of Altarelli—Parisi’s P-function. 1(0) stands for NLL(LL). NLAL: Selection of the form of the strong coupling constant, a 8. When it is 1(0), a8 is calculated in the form of NLL (LL). 2 for KZAP: =Selection the for argument strong of partons. It is Qto the KZAP 0, x(l of x)Q2 KZAP of = the 1, and (1 coupling x)Q2 forconstant KZAP =in2.the Theevolution last two cases correspond squared transverse momentum. KZQQG: Selection of the argument of the strong coupling constant at the primary q4g-vertex. The meaning is the same as for K.ZAP. —

—

90

K. Kato, I Munehisa

/

NLLjet: a MC code fore ~e - QCD jets

NSHUFI, NSHUF2: (Unused now.) NLMODE: When it is 1(0), the evolution is done in NLL (LL) mode. NAO: When it is 1(0), the condition of angular ordering is imposed (not imposed). /TPNEGG/ PNEGG(O:NIDIM), DELTG, NTMXG, NINTG Table of non-branching probability, HNE for gluon. PNEGG: PNEGG(K) is log HNE(Q~SX, Q~) for gluon. Q~axis the squared center-of-mass energy and logQ~= log Q~nm+ K x DELTG. The values of this table are negative and PNEGG(NTMXG) = 0 by 2) log(Q~ 2) are found in /LIMKIN/ as QML and definition. It should be noted that log(Q~~~/Aand 1~/A QOQGL. DELTG: log(Q~~ 8/Q~j~)/NTMXG, where Q~1is the available minimum virtuality. NTMXG: The size of the table. NINTG: (Unused). /TPNEQQ/ PNEQQ(O:N1DIM, 10), DELTQ(1O), NTMXQ, NINTQ The same tables as in /TPNEGG/ but for quarks. The second index of PNEQQ and the index of DELTQ specify the flavor of quark. /TPNEGQ/ PNEGQ(0:N1DIM), DELGQ, NTMGQ, NINGQ The same tables as in /TPNEGG/ but for the branching of a gluon 2 into In order to keep are q4. identical to those in consistency, the division of the table, the maximum and minimum of Q /TPNEGG/. /PRVTX/ S2MOD(1O), FMODE(10) Fractions for selecting the quark flavor and the primary vertex. S2MOD: The fraction of the two-jet, that depends on the quark flavor. FMODE: The flavor fraction. /COLOR/ CF, CA, TR, PT, BETAO, BETA1 Values for color factors and ir. See eq. (A.ll). CF: CF. CA: CA. TR: TR. P1: ¶=3.14l59... BETAO: The lowest-order coefficient of the beta function, /3~= ~-(IIc/A of flavors. BETA1: The next-to-leading-order coefficient of the beta function, ~

—

=

4TRNf), where N~ is the number

~(34C~

—

(2OCA + I2CF)TRNf).

/LIMKIN/ QM(O:l0), QOM(O:1O), QO2M(O:1O), QOQG2(O:1O), QOQQ2(IO), QOQGL(0:1O), QML Values for the minimum virtual mass and related quantities. Index of arrays is 0 for gluon, and 1, 2, 3, 4, 5 and 6 for the u-, d-, s-, c-, b- and t-quark, respectively. QM: Mass of quark, in GeV. QM(0) = 0. QOM: The lowest virtual mass, ~/Q 2, in GeV. QOM2: Square of QOM, in GeV2. 02 + QM QOQG2: The minimum virtual mass squared for a system of a quark and a gluon, in GeV2.

K. Kato, T.

Munehisa /

91

NLLjet: a MC code fore ~e — QCDjets

QOQQ2: The minimum virtual mass squared for a system of a quark and an antiquark, in GeV2. QOQGL: log QOQG2/A2. QML: log s/A2 where s is the center-of-mass energy squared. /CWEAKO/ NWEAK NWEAK: It determines the quark-current coupling. 1 (0) stands for the electro-weak theory (QED). /CWEAK2/ ZMASS, WMASS, HMASS, SIN2W, ZGAMM, WGAMM, WKCP2, COEF1X(2), COEF2X(2), COEF3X(2), COEF4X(2), COEF5X(2), COEF6X(2) Physical parameters for electro-weak theory. Some of these are set in QCINPT and QCDATA as input values and the remainders are set in TBFMD. At present, WMASS, HMASS and WGAMM are not used in NLLjet. ZMASS, WMASS, HMASS: Masses of Z°,W and Higgs particle. SIN2W: x~,=sin29~,. ZGAMM, WGAMM: Width of Z°and W. WKCP2: 1/16 sin29~cos29~. COEFxX(J): It represents coefficients for two-jet distribution. J = 1 represents the up-state of weak isospin and J = 2 the down-state. These are used in the package of radiative correction. /DSTRJT/ RQMS(iO), RGMS(1O), COEF(6, 10), F2J(4, 6), P1 Coefficients for angular distributions of three-jet and two-jet. The index of the above arrays represents the flavor of quark. RQMS(I): 2mq/%/~for the I-th quark. RGMS: 2Q 0/V~for all quarks. 9/i6IT)v(Cvvfvv + COEF: The differential crossofsectionwhere of two-jets for the I-th quark is da/d~7 ( CAAfvv + CvAfvA) in units C~= COEF(1, I), CAA = COEF(2, I),= CvA = COEF(3, I) + COEF(4, I), fvv = 2 + v2 cos20, fAA = v2 + v2 cos20, fVA = 2v cos 9, and v = ~/i 4m~/s. F2J: The differential cross section of two-jets for the I-th quark is given by the function FJET2. FJET2 = FJ2(l, I) + FJ2(2, I) cos 0 + FJ2(3, I) cos20 and the maximum of FJET2 is FJ2(4, I). P1: IT=3.l4i59... —

—

/QCRECV/ NRECOV NRECOV: It is referred in QCPREP and controls the save/load action of variables in /TPNEGG/, /TPNEQQ/, /TPNEGQ/ and /PRVTX/. When it is 0, the data file is not used. When it is negative, the tables are written on a file of unit number= —NRECOV. When it is positive, the tables are read from a file of unit number = NRECOV. /QCRND/ NSEED NSEED: Value for seed of the random number. /COMAL1/ CO, Ci, FL2, FL2M, QCUT CO: 13 4Tr. Cl: $ 0/4IT$ 1/ 2. 0. FL2: A A2 that is defined as a~(Q~) FL2M: = QCUT: Q~.

2).

4IT//3 0

log(Q~/A

K. Kato, T. Munehisa

92

/ NLLjet:

a MC code fore e - QCDjets

/COMAL2/ ff1, B2, B3, B4, BS These values are the coefficients of asymptotic form of the next-to-leading quark P-function. At

present, b1(= Bi) is used as Pt1~(x) b1/(l —

—

x) as x

—~

1.

/COMPES/ NSAME(10) NSAME: If mass of the I-th quark and that of J-th quark (J> I) are the same, we set NSAME(J) = I to skip the calculation of the Sudakov factor for the J-th quark. Otherwise, NSAME(J) = 0. /COMPEW/ WG1(0:Nl), WG2(0:N2), WG3(0, N3) Weights for numerical integration by the Simpson method. /ALPHAL/ QCUT, ALPHMX To set the upper limit of the strong coupling constant. 2 for the argument of the strong coupling constant. QCUT: The The minimum Q value of the strong coupling constant. ALPHMX: maximum /CMSANG/ CO(3), SI(3) These keep the direction of primary partons. They are set in CMSPRV and used in CMSROT to rotate the partons’ four-momenta. /CKIN1/ PIP2, P1P3, P2P3, El, E2, E3, P1, P2, P3, C12, C13, C23 These values give p, p1 j, E,, I i-’ I~ and cos 0,. in q4g primary vertex. They are set in CHKKN1.

/RCCOM1/, /RCCOM2/, /RCCOM3/, /RCCOM4/, /RCCOM5/ These five common areas are used by the package for radiative correction. /SV3GWZ/ Z, ZI This common area is used by SV3GW to transfer variables to some low-level routines. /QCWARN/ NWARN NWARN: (Unused now.)

Appendix D. Input data NLLjet is controlled by a set of input parameters. In this appendix, we describe how to specify them. Input data are read from unit FTO5FOOI and consist of a number of lines. Each line consists of a keyword and of a numerical parameter or a string parameter. For some keywords, there is no parameter. The number of parameters is one, except for the QMASS line. The keyword and parameters are to be separated by spaces. The location of a keyword and parameters in a line is arbitrary. The ordering of keywords is also arbitrary, except for the QMASS line. Since each parameter has its default value, only values that the user wants to change are to be given. Below we briefly present keywords and their meaning. The keywords must be types in capitals as shown. After a keyword, we specify the parameter type by [I], [R], [C], or [X] to denote that the type is integer, floating, or string, or that there is no parameter for the keyword. The value quoted in parenthesis is the default value. TITLE [C] The title string. (String to specify date) ROOTS [R] Center-of-mass energy, in GeV. (50)

K Kato, T. Munehisa

Q02 [R] LAMBDA [R] NFLAV [I] QMASS [R]

NRECOV [I]

/ NLLjet: a

MC code fore ‘e — QCDjets

93

Minimum virtuality of partons, in GeV2. (1.0) A~,in GeV. (0.2) Number of quark flavors. (5) Quark masses, in GeV. This keyword is exceptional for two reasons. First, it must appear after the NFLAV line if this line exists. Second, the number of parameters is not 1 but NFLAV. Ordering of quarks is u, d, s, c, b and t. (0.0, 0.0, 0.0, 1.5, 5.0, 30.0) Specify the unit number of the file for storing the Sudakov-factor table and the fraction of jets. This is to save computer time in the initialization stage. If it is 0, then no read/write action is done for the file. If it is negative, then calculated values are written on the file of unit number = NRECOV. If it is positive, the step to calculate tables is skipped and the values are read from the file of unit number = NRECOV. (0) Number of events to generate. Please note that NEVEN and the following NPRINT are used in the standard main program that is distributed with NLLjet. When the user employs his/her own main program, these parameters are meaningless. (10) Interval for printing the list of generated events. Events are printed once in NPRINT event. (10) Parameter for separating primary two- and three-jets. (0.1) To order that the quark current only couples to the photon and not to Z°.This is not the default status. Mass of Z°,in GeV. (91.1) Width of Z°,in GeV. (2.5) sin20~.(0.223) Seed number for random number generation. When the user modifies the routine for random numbers to suit the host machine, this value is meaningless. (1234567879) To order that the leading-logarithmic mode is used. This is not the default status. To order that the application of the angular-ordering condition is suppressed. This is not the default status. —

NEVEN [I]

NPRINT [I] DELTA [R] QED [X] ZMASS [R] ZGAMM [R] SIN2W [R] NSEED [I] LL [X] NOAO [X]

The following parameters control the radiative correction. IRC [I] Flag to switch radiative correction. When IRC = 0, radiative correction is not taken into account. When IRC = 1 or 2, correction for initial radiation is included. The difference between IRC = 2 and IRC = 1 is that the former uses the formula in which soft photons are exponentiated. (0) EPMIN ER] Minimum of photon energy, in GeV. No meaning for IRC = 0 or IRC = 2. (1.0) EPMAX ER] Maximum of photon energy, in GeV. No meaning for IRC = 0. (ROOTS* 0.95/2.0) BSPREAD [I] Flag for the width of beam energy, in GeV. BSPREAD = 0/1 is off/on. (0) SIGMO [RI Width of beam energy, in GeV. (0.0)

Appendix E. Example of the main program and JCL for the FACOM at KEK The following programs are two examples of a main program for generating events by NLLjet. After calling QCINPT and QCPREP, QCGEN is called repeatedly. C C C C

Sample main program for NLLjet (1) Generate events NEVEN times and print them once for NPRNT events.

94

K. Kato, I Munehisa

C C C C C C

IMPLICIT REAL ...

...

*

*

8(A-H, O-Z)

Read parameters CALL QCINPT (WBEAM, NEVEN, NRADC, NPRNT) Prologue for NLLjet CALL QCPREP

If you want to do some initialization, do it here. ...

C C

C C C C C C C C

a MC code for e e — QCD jets

(default: NEVEN = 10, NPRNT = 10) No radiative correction, No hadronization.

C C C C C C

/ NLLjet:

...

Event Loop DO 100 1=1, NEVEN Generate one event IFLAV =0 CALL QCGEN (IFLAV)

*

If you want to hadronize partons, call hadronization routines here.

*

Now one event generated. If you want to analyse it, do it here.

*

*

If you do not need to print, delete following lines. Print out parton list NEVE =1 NCHK = 1 IF(MOD(I,NPRNT).EQ.O) CALL QCOUT (NEVE, NCHK, WBEAM) ...

C 100 CONTINUE C C C

*

NEVEN events generated. If you want to do summary, do it here. STOP END

If the user wants to generate jets of a specified flavor, it can be done by setting IFLAV in the main program. When the user wants to use the radiative correction package in NLLjet, the main program is as follows. We need the common area /NLLPAR/ and PHOT(4) for the photon’s momentum. QCRCGN generates the photon for a radiative event and it also determines the flavor IFLAV. Since the development of the parton shower is done in the center-of-mass system of reduced energy, the system of partons is transformed into the laboratory system by QCRCBS. We also shift 6 for a massive quark since it is impossible to generate primary three-jets of massive quark when the radiated photon takes away a large amount of energy.

K Kato, T. Munehisa

C C C C C C C C

/

NLLjet: a MC code fore ~e — QCDjets

Sample main program for NLLjet (2) Generate events NEVEN times and print them once for NPRNT events. (default: NEVEN = 10, NPRNT = 10) With/without radiative correction.

IMPLICIT REAL C C C C

*

8(A-H, O-Z)

PHOT(4) = 4-momentum of radiative photon IPHOT = 0... Non-radiative event = 1 Event with a radiative photon DIMENSION PHOT(4) COMMON /NLLPAR/ ROOTSO, ROOTS, Q2MAX, Q02, FLAMD, FLMD2 DELTA(1O), QMASS(l0), SCHEME, NFLAV ...

*,

C C

...

C C

...

C C C C C C

*

Read parameters CALL QCINPT (WBEAM, NEVEN, NRADCO, NPRNT) Prologue for NLLjet CALL QCPREP

If you want to do some initialization, do it here. ...

...

Event Loop DO 100 I =1, NEVEN Radiative correction, generate photon IF (NRADCO.NE.0) THEN CALL QCRCGN (ROOTSX, ROOTS, PHOT, IPHOT, IFLAV)

C IF (IPHOT.NE.O) THEN DELSAV = DELTA(IFLAV) DELTA(IFLAV) = MAX(DELTA(IFLAV), 2 END IF *

*

SQRT (QMASS(IFLAV)

C END IF C C

...

C C

...

generate partons CALL QCGEN (IFLAV) boost system if initial state photon is radiated. IF (NRADCO.NE.0) THEN CALL QCRCBS (ROOTSX, PHOT, IPHOT) IF (IPHOT.NE.0) DELTA(IFLAV) = DELSAV END IF

*

*

+ Q02)/ROOTS)

95

96

C C C C C C C C

K. Kazo, T. Munehisa

/

NLLjet: a MC code fore + e — QCDjets

*

If you want to hadronize partons, call hadronization routines here.

*

Now one event generated. If you want to analyse it, do it here.

*

If you do not need to print, delete following lines. Print out parton list and hadron list IF (MOD(I,NPRNT).EQ.O) THEN NEVE = 1 NCHK = 1 CALL QCOUT (NEVE, NCHK, WBEAM) END IF

*

...

C 100 CONTINUE C C C

NEVEN events generated. If you want to do summary, do it here.

*

STOP END Below we show an example of job-control statements to run NLLjet on the FACOM at KEK. Here XXXX is a user-id and XXXX.NLLJET.FORT is a partitioned dataset that stores the FORTRAN source of NLLjet. The member NLLJET consists of 99 modules described in appendix B, and MAINO is the first main program given above. As a stream file, FTO5FOO1 is defined to specify control parameters that are explained in appendix D. Since NRECOV is defined as a control parameter, the file [email protected] is declared in order to save the Sudakov-factor table and that of the fraction of the primary vertex. //XXXX JOB CLASS = M //

*

NLLjet: JOBO: Generate Parton Shower //* //EXEC FORT7CLG, // PARM.FORT = ‘OPT(3), SOURCE’, // PARM.LKED = ‘NOMAP, NOLIST,’ //FORT.SYSIN DD DSN = XXXX.NLLJET.FORT(MAINO), DISP = SHR // DD DSN = XXXX.NLLJET.FORT(NLLJET), DISP = SHR //GO.FTO5FOO1 DD TITLE TEST RUN ROOTS 29.0 LAMBDA 0.20 DELTA 0.1 Q02 1.0 NFLAV 5 QMASS 0.325 0.3250 0.5 1.6 5.0 NEVEN 10000 NPRINT 1000 //

*

*

K. Kato, I Munehisa

/

NLLjet: a MC code fore ~e — QCDjets

97

NRECOV -30 IRC 0 /*

//GO.FT3OFOO1 DD DSN = [email protected],DISP = (NEW,CATLG), // SPACE = (TRK,(2,2),RLSE), DCB = (RECFM = VBS,BLKSIZE = 23689) //

References [1]T. Muta, Foundations of Quantum Chromodynamics (World Scientific, Singapore, 1987). K. Hagiwara et al., Suppi. Progr. Theor. Phys. 77 (1983) 1. [2]B. Bambah et al., preprint CERN-TH 5466/89, in: Z Physics at LEP 1, vol. 3, ed. G. Altarelli, R. Kieiss and C. Verzegnassi, CERN 89-08, Geneva (1989) p. 143. [3]M. Bacé, Phys. Lett. B 78 (1978) 132. [4] K. Komshi, A. Ukawa and G. Veneziano, Nucl. Phys. B 157 (1979) 45. [5] R. Odorico, Nucl. Phys. B 172 (1980) 157. G.C. Fox and S. Wolfram, Nucl. Phys. B 168 (1980) 285. T.D. Gottschalk, NucI. Phys. B 214 (1983) 201; 227 (1983) 413. [6] G. Marchesini and B.R. Webber, Nuci. Phys. B 238 (1984) 1. B.R. Webber, Nucl. Phys. B 238 (1984) 492. [7] T. Sjostrand, Comput. Phys. Commun. 39 (1986) 347. T. Sjostrand and M. Bengtsson, Comput. Phys. Commun. 43 (1987) 367. [8] J. Kalinowski, K. Konishi and T.R. Taylor, Nucl. Phys. B 181 (1981) 221. J. Kalinowski, K, Konishi, P.N. Scharbach and T.R. Taylor, Nucl. Phys. B 181 (1981) 253. [9] J.F. Gunion and J. Kalinowski, Phys. Rev. D 29 (1984) 1545. J.F. Gumon, J. Kalinowski and L. Szymanowski, Phys. Rev. D 32 (1985) 2303. [10] K. Kato and T. Munehisa, Mod. Phys. Lett. A 1 (1986) 345; Phys. Rev. D 36 (1987) 61. [11] K. Kato and T. Munehisa, Phys. Rev. D 39 (1989) 156. [12] V. Sudakov, Soy. Phys. JETP 3 (1956) 65. [13] G. Altarelli and G. Parisi, Nuci. Phys. B 126 (1977) 298. [14] A. Bassetto, M. Ciafaloni and G. Marchesini, Phys. Rep. 100 (1983) 201 and references therein. [15] T.D. Gottschalk and M.P. Shatz, Caltech preprint, CALT-68-1172; Phys. Lett. 150 B (1985) 451. [16] T. Kamae, in: Proc. 24th mt. Conf. on High-Energy Physics, Munich, 1988, eds. R. Kotthaus and J.H. Kuhn (Springer, Berlin, 1989) p. 156. [17]K. Kato and T. Munehisa, KEK Report 84-18 (1984); 87-5 (1987). [18]K. Abe et al., VENUS Collaboration, Phys. Lett. B 240 (1990) 232. [19]F.A. Berends, R. Kleiss and S. Jadach, Nucl. Phys. B 202 (1982) 63.