FTUV110130 KATP102011 LPN1120 SFBCPP1122
Towards at NLO QCD: Bosonic
[2mm] contributions to triple vector boson
[2mm] production plus jet
[2mm]
Francisco Campanario
Institut für Theoretische Physik,
Universität Karlsruhe, KIT,
76128 Karlsruhe, Germany
In this work, some of the NLO QCD corrections for are presented. A program in Mathematica based on the structure of FeynCalc which automatically simplifies a set of amplitudes up to the hexagon level of rank 5 has been created for this purpose. We focus on two different topologies. The first involves all the virtual contributions needed for quadruple electroweak vector boson production, i.e. . In the second, the remaining “bosonic” corrections to electroweak triple vector boson production with an additional jet () are computed. We show the factorization formula of the infrared divergences of the bosonic contributions for VVVV and VVVj production with . Stability issues associated with the evaluation of the hexagons up to rank 5 are studied. The CPU time of the FORTRAN subroutines rounds the 2 milliseconds and seems to be competitive with other more sophisticated methods. Additionally, in Appendix A the master equations to obtain the tensor coefficients up to the hexagon level in the external momenta convention are presented including the ones needed for small Gram determinants.
August 8, 2021 —19:46
1 Introduction
The incoming data from the CERN Large Hadron Collider (LHC) will require precise theoretical predictions for a variety of signal and background processes. Processes with two vector bosons are of vital importance since they constitute background signals to Higgs and top physics as well as to physics Beyond the Standard Model (BSM) (for an overview see e.g. Ref. [1]). A tremendous effort has been carried out by the scientific community to compute these processes accurately. NexttoLeading Order (NLO) QCD diboson production was computed in Refs. [2, 3, 4], and with an additional jet in Refs. [5, 6, 7, 8, 9, 10, 11]. Electroweak (EW) diboson production in association with two jets was computed in Refs. [12, 13, 14]. Recently, the NLO QCD corrections to [15] and [16] production were computed using the generalized Ddimensional unitarity framework for the calculation of the one loop virtual amplitudes. An impressive progress in the calculation of multiparton oneloop amplitudes has been achieved not only applying new techniques [17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28, 29, 30, 31, 32, 33, 34, 35, 37, 38, 39, 40, 41, 42, 43, 44, 45, 46, 36, 47], but also based on traditional methods [48, 49, 50, 51, 52, 53, 54, 55, 56, 57, 58, 59, 60, 61, 62, 63, 64, 65, 66, 67, 68, 69, 70, 71, 72, 73, 74]. In this paper, relying on traditional techniques, we compute some of the virtual corrections needed for production. We focus on some of the topologies that appear, namely QCD “bosonic” oneloop corrections to the diagrams contributing to and production. The strategy followed is to collect Feynman diagrams with a given topology in groups which can be easily checked and reused in other processes. Thus, our aim is not only to provide results for production, including a second calculation of production using a different method, but also to provide a set of routines that can be used for many other interesting processes, such as at NLO QCD [75]. It is known that in the calculation of oneloop multileg amplitudes, the presence of small Gram and Cayley determinants might yield unstable results. In this paper, we show that the use of higher precision in the numerical determination of the tensor integrals (similarly to Ref. [76]) together with the use of Ward Identities to identify the unstable points, and to a minor extent, following Ref. [57], the use of special routines for the determination of small Gram determinants for the and functions, solve this problem.
The paper is organized as follows; in Section 2, the method used to perform the calculation, the different contributions computed and the tests performed to guarantee the correctness of the results are described. In Section 3, a discussion of the instabilities and the timing of the amplitudes computed are shown. The conclusions are given in Section 4. In Appendix A, the tensor reduction routines used are presented. In Appendix B, numbers for the contributions involving the hexagons are given, including proofs of the factorization of the infrared divergences and other additional numerical tests. Finally, in Appendix C, we give the color factors needed for testing the amplitudes.
2 Calculational details
In the calculation of the NLO QCD virtual corrections to , different subprocesses contribute; processes involving two quark pairs, e.g., and processes with one quark pair, e.g., ^{1}^{1}1Note that the NLO QCD corrections for production involve only subprocesses with two quark pairs.. The latter can be separated into “one loop fermion” corrections and “bosonic” corrections, Fig. 1.
Among the “bosonic” contributions, different topologies appear, Fig. 2. We will focus on diagrams with the first two topologies, i.e., one loop QEDlike corrections and diagrams formed with one triple gluon vertex. The remaining ones will be discussed in future publications.
To compute these amplitudes, a function in Mathematica [77] using FeynCalc [78] which automatically simplifies a set of amplitudes up to Hexagons of rank 5 has been created. Throughout the calculation, the quarks are considered to be massless, the anticommuting prescription for is used and we work in the Feynman Gauge. The result is given in terms of tensor integrals following the PassarinoVeltman convention [48], Appendix A, and written automatically in to FORTRAN routines. Nevertheless, the amplitudes can be evaluated also in Mathematica with unlimited precision which is used for testing purposes. To achieve that, the scalar integrals, the tensor reduction formalism to extract the tensor coefficient integrals, and also the helicity method described in Refs. [79, 80] to compute the spinor products describing the quark lines of Figure 2 have been implemented at the FORTRAN and Mathematica level. For the determination of the tensor integrals up to the box level, we have implemented the PassarinoVeltman tensor reduction formalism [48]. We have also applied the LU decomposition method to avoid the explicit calculation of inverse Gram matrices by solving a system of linear equations which is a more stable procedure close to singular points. Finally, for singular Gram determinants, special tensor reduction routines following Ref. [63] have been implemented, however, the external momenta convention (Passarinolike) was used. The impact of these methods is discussed in detail in Section 3. For pentagons, in addition to the PassarinoVeltman formalism, the method proposed by Denner and Dittmaier [63], applied also to hexagons, has been implemented. For that, the recursion relations of Ref. [63] in terms of the PassarinoVeltman external momenta convention have been rederived. This last method is used for the numerical implementation at the FORTRAN level in Section 3. All the recursion relations used can be found in Appendix A.
In the following, a more detailed description of the method used is given. To extract the rational terms, two simplifications are done at the Mathematica level. To illustrate them, the amplitude of a simple vertex diagram is used,
First, a simple Dirac reordering manipulation is applied to explicitly contract
repeated indices and to obtain all the terms proportional to the space dimension, , coming from
contractions, i.e.,
(2.1) 
Second, all Dirac’s structure containing the loop momenta is pulled to the right, such that resulting terms like are canceled against one of the integral denominators,
(2.2)  
where is the scalar two point function defined correspondingly to Eq. (A.1). In this way we avoid possible additional terms proportional to coming from the Lorentz structure of the tensor integral. More concretely, the same term could result in terms proportional to :
(2.3)  
with , a three point tensor coefficient integral of the tensor integral , defined by Eq. (A.6), and obtained recursively from Eq. (A). Further simplifications are obtained using the Dirac equation of motion and rewriting the pair as a difference of two propagators, e.g., , which are canceled against the denominators. The remaining rational terms stem from ultraviolet divergent tensor coefficients, which are treated independently within the tensor reduction routines. As an example, the finite contribution of the ultraviolet divergent tensor coefficient , for massless propagators, following Eq. (A) is obtained by,
(2.4) 
where we have taken and series expand in all the scalar and tensor coefficient integrals, e.g.,
(2.5) 
with . After these steps, the amplitudes, , in terms of scalar and tensor coefficient integrals can be written by,
(2.6) 
where is the amplitude that one would obtain performing the Dirac algebra manipulation in four dimensions, , and contains the rational terms and vanish in Dimensional Reduction (). To get this expansion, we have only rewritten the space dimension D as such that for example Eq. (2.2) is given by,
(2.7) 
Both and are decomposed in terms of,
(2.8) 
where SM is a basis of Standard Matrix elements corresponding to spinor products describing the quark line of Figure 2 which are computed following the helicity method [79, 80] with a defined helicity, . F1 are complex functions which are further decomposed into dependent and independent loop integral parts,
(2.9) 
is a monomial function at most for each polarization vector . contains kinematic variables (), the scalar integrals ()^{2}^{2}2Note that throughout the paper to name the different npoint function integrals, , the alphabeticallyordered labeling frequently found in the literature is also used., and the tensor integral coefficients . The latter obtained numerically from the scalar integral basis following the recursion relations of Appendix A. To illustrate this notation, the following example is given,
(2.10) 
is the positive helicity projector, . We note here that this parametrization is convenient for example to evaluate the amplitude for different polarization vectors or for performing gauge tests since some of the functions will remain unchanged. Finally, the full result is obtained from and using the finite and the coefficients of the poles of the scalar and tensor coefficient integrals (see e.g. Eq. (2.5)). From , we get
(2.11) 
where is the finite contribution obtained using the finite pieces of the scalar and tensor coefficient integrals including the finite contributions from rational terms arising in ultraviolet tensor coefficient integrals, Eqs. (2.4) and (A). and are obtained from the and pole contributions, respectively. Similarly, from , one gets,
(2.12) 
where and are obtained from the and poles, respectively. Indeed, is analytically zero since rational terms in oneloop QCD amplitudes, omitting wave renormalization graphs (WRF) are of UV origin [81]. Eq. (2.12) will allow us to check this statement numerically. Additionally, up to the pentagon level, the divergent part is computed in Dimensional Regularization analytically in Mathematica. A library containing all the divergent tensor coefficients for massless particles including pentagons of rank 4 has been created with this purpose. This allows us to obtain also the terms of Eq. (2.12) and the divergent parts of Eq. (2.11) analytically up to this level.
It is known that the IR divergences depend on the kinematics of the external particles involved in the process, i.e., whether they are massless/massive onshell/offshell particles. Nevertheless, since we reconstruct the divergences numerically, to compute generally the amplitude of a particular diagram for a given helicity, , it is sufficient to consider offshell vector bosons, in Figure 2. This basically means that all terms are kept and the transversality property of the vector onshell bosons () is not applied analytically, allowing us, the latter, to perform gauge tests () for onshell massive particles. We note here that the polarization vectors can be understood as generic effective currents, , which can include the leptonic decay of the vector bosons or physics BSM. Thus, denoted by the amplitude of the first diagram of Figure 2 would be given by,
(2.13) 
where is the strong unrenormalized coupling, is the corresponding diagram with offshell vector bosons (or effective currents) with color indices and given in terms of Eq. (2.6) which is computed only once and for all vector bosons, . is a color diagram dependent factor, e.g., and from now on the color subindices will be omitted. is the coupling following the notation of Ref. [80] with e.g., , , , where is the weak mixing angle, is the third component of the isospin of the (lefthanded) fermions, and the electric charge.
This procedure is advantageous for several reasons. First, when considering offshell vector bosons, for example, diagrams with the same QEDlike topology as the first one of Fig. 2 will give us, already, all the corrections needed for EW quadruple vector boson production, and part of the and contributions, independently whether we are considering massless or massive particles or the bosons emitted from the quark line are gluons. will associate the correct color factor to each diagram depending whether one is considering , or and the corresponding couplings. Also, the leptonic decay of the vector bosons and new physic effects can be immediately incorporated by using the appropriate effective currents. Second, cross term related diagrams are obtained from the same analytical amplitude just by permuting the vector boson momenta and polarization, and thus, e.g., from the 24 cross related diagrams to the first one of Fig. 2, only the permutation depicted has to be computed. Moreover, additional checks can be implemented, as for example the factorization of the IR divergences against the born amplitude for different processes starting from the same analytical structure. In the following, we will refer to the first topology of Fig. 2 as contributions to and, to the second as “bosonic” contributions to since they appear for the first time in these processes.
2.1 Contributions to
In this section, all the “bosonic” QEDlike loop corrections along a quark line needed for production are considered. These contributions can be classified by I) virtual corrections along a quark line with two vector bosons attached (first diagram of Figure 3) II) virtual corrections with three vector bosons attached to the quark line (second and third diagrams of Figure 3) and III) virtual corrections along a quark line with four vector bosons attached (last diagram of Figure 3).
The strategy followed is to search for a minimal set of universal building blocks from which every amplitude can be constructed. We use the effective current approach described and applied in Refs. [7, 10, 11, 82, 83, 84, 85, 86]. As illustration, the first diagram of Fig. 3 is used. This can be written as,
(2.14) 
where the color indices have been omitted. Here, and represent effective polarization vectors including finite width effects in the scheme of Refs. [87, 12] and propagator factors, e.g.,
(2.15) 
with , the width of the vector boson, and , the triple vertex which can also contain the leptonic decay of the EW vector bosons including all offshell effects or BSM physics. In this manner, we can treat diagrams with triple vertexes as external offshell legs (in the sense that they will have the same IR behavior). We can then concentrate in computing, instead of , the virtual correction to two massive vector bosons attached to the quark line, , or equivalently , where the polarization vectors or effective currents have been factored out. Thus, the strategy will be to combine in groups all the virtual corrections to a given bornamplitude configuration, independently of the effective currents or polarization vectors attached to the quark line. Furthermore, the order of the gauge bosons are fixed and the full amplitude will be recovered by summing over the physicallyallowed permutations. Thus, three universal virtual contributions are left, corrections to a born amplitude with two, three and four vector bosons attached to the quark line for a given order permutation.
I) The virtual corrections to the Feynman graph with two vector bosons and attached to the quark line for a given orderpermutation with incoming momenta and are depicted in Figure 4
with kinematics,
(2.16) 
The sum of the four amplitudes with a given helicity is written in terms of,
(2.17) 
where is called from now on “boxline” contribution. Although, we are interested on offshell vector bosons emitted from the quark line (effective polarization vectors), the divergent contributions and were also computed analytically for the different kinematic configurations. These configurations are obtained by considering the vector bosons massless or massive and all possible ordering permutations, i.e., . From the analytical calculation, it is confirmed that the rational terms arise from diagrams with ultraviolet divergences, thus, is zero. For configurations with massless vector bosons, , the transversality property of the bosons must be used, . A general proof of this statement was presented in Ref. [81]. Moreover, is zero since there are not UV poles. For the selfenergy diagram, =0 since the divergences are only . Finally, for , there are not collinear and soft singularities simultaneously given rise to poles, so that, =0.
II) The virtual corrections to the Feynman graph with three vector bosons , and attached to the quark line for a given permutation are depicted in Figure 5 with kinematics,
(2.18) 
The sum of the eight diagrams with a given helicity is written in terms of,
(2.19) 
where is called from now on “penline” contribution. We have computed also the divergent contributions and analytically for the eight different kinematic configurations, . From the analytical calculation, it is verified that once the transversality property for massless onshell particles is applied, the rational terms arise from diagrams with ultraviolet divergences. Therefore, and are zero. Analogously to the boxline, for selfenergies, =0 as well as =0.
III) The virtual corrections to the Feynman graph with four vector bosons , , and attached to the quark line for a given permutation are depicted in Figure 6 with kinematics,
(2.20) 
The sum of the thirteen diagrams with a given helicity is written in terms of,
(2.21) 
will be called “hexline” contribution in the following. For the hexline contribution, the divergences are not computed analytically. The explicit and naive recursive reduction and simplification of the tensor integrals for hexagons in terms of poles exceeds the capacity of Mathematica with 4 GB of memory RAM. Therefore, Eqs. (2.11,2.12) are used to obtain numerically the different contributions both at the Mathematica and FORTRAN level. We have checked with up to 30000 digits of precision in Mathematica for the 16 different kinematic configurations, i.e., and several phase space points that is zero for all possible kinematic configurations (for massless particles the transversality property must be used) as well as . At the FORTRAN level for nonsingular points, the proof works at the working precision level.
2.1.1 Checks
In this section, the tests performed in order to ensure the correctness of the calculation of the 1 loop diagrams/contributions are explained. We consider electroweak vector boson production (, (W,Z,). Then, the color factor for all diagrams and contributions is proportional to . This is enough to test not only the individual diagrams but also the different contributions. For electroweak vector bosons, the divergent terms for the boxline and penline are known analytically, the general form generalizes for n vector boson emission and reads,
where is the corresponding Born amplitude with helicity and the square of the partonic centerofmass energy. If we use in Eq. (A.2) for the definition of the scalar and tensor integrals and , then, the prefactor of Eq. (2.1.1) is reproduced and each of the terms of Eq. (2.1.1) can be identified, for example, for the hexline contribution, Eq. (2.1), with the different finite and pole terms as,
(2.23) 
These relations are very important for testing purposes. Note that the left hand quantities for each of the above lines are obtained numerically from the same analytical expression, Eqs. (2.11,2.12). They only differ in which terms of the expansion of the scalar and tensor integrals are used. Therefore, the numerical check of the factorization of the singularities and rational terms provides a strong check of the correctness of the finite terms . We have checked the factorization formula analytically for the boxline and penline contributions for all different kinematic configurations. To cast the singularities in this form for massless onshell vector bosons, (), the transversality property of the particle must be used. Otherwise, additional singularities proportional to the product appear, with the onshell massless particle momentum. Thus, we have checked that the presence of additional massless particles does not introduce new singularities or new rational terms in agreement with Ref. [81]. For the hexline contribution, we have checked numerically with high precision in Mathematica the factorization formula for the 16 different kinematic configurations and for different phase space points . The coefficients multiplying the poles of Eq. (2.1.1) are obtained with up to 30000 digits of precision at least (at the FORTRAN level, for nonsingular points, the proof works at the working precision level). Nevertheless, we will assume the result to be analytic due to the high precision achieved. An analytical proof can be obtained using the method described in Ref. [88] which is not automatized within our approach.
Concerning the origin of the rational terms in Eq. (2.1.1), we note that for electroweak boson production, the wave renormalization functions (WRF), which are zero in Dimensional Regularization for massless onshell quarks, together with are UV finite. Therefore, all the divergences in Eq. (2.1.1), after adding the WRF, become of IR origin including the terms containing the rational factor, , of UV origin. This is clear since the WRF are zero due to the cancellation of IR and UV poles. Treating the UV and IR divergences separately, one observes that the UV poles of the WRF cancel the rational terms of Eq. (2.1.1) and only IR divergences remain (see also Ref. [81]).
The factorization proof already provides an important check of the correctness of the calculation. We can use the factorization formula to perform an additional test. The factorization of from the scalar integrals introduces an additional dependence in on this variable. If we consider (s) as an independent scale energy variable, , and series expand in in the second piece of Eq. (2.1.1), we obtain a new finite term, , of the form ^{3}^{3}3Note that and , which is , have the same analytical structure. They only differ in the input scalar integrals.,
(2.24) 
such that the terms compensate the dependence propagated in all the scalar and tensor coefficient integrals in the complex . Thus, is independent. This fact was satisfied for nonsingular points in the three contributions at the working precision level in the FORTRAN code (for the hexline see Tab. 12) and at least with 30000 digits in the Mathematica code. Moreover, we have implemented Ward identity tests for the virtual corrections in FORTRAN and Mathematica at different levels of complexity:

At the level of single diagrams.

For the , , contributions.

At the level of gauge invariant quantities.

Subset of amplitudes invariant for a specific replacement ().

Specific contractions that make the contributions to vanish.
1) Under the replacement of the polarization vector by its momentum (), relations among the different diagrams/topologies can be found. As illustration, one can consider the hexagon of Fig. 7,
(2.25) 
where . Contracting one of the open indices by the corresponding momentum and expressing the contracted gamma matrix as the difference of two adjacent fermionic propagators, the hexagon is reduced to a difference of two pentagon integrals,