Automated one-loop calculations with GoSam

We present the program package GoSam which is designed for the automated calculation of one-loop amplitudes for multi-particle processes in renormalisable quantum field theories. The amplitudes, which are generated in terms of Feynman diagrams, can be reduced using either D-dimensional integrand-level decomposition or tensor reduction. GoSam can be used to calculate one-loop QCD and/or electroweak corrections to Standard Model processes and offers the flexibility to link model files for theories Beyond the Standard Model. A standard interface to programs calculating real radiation is also implemented. We demonstrate the flexibility of the program by presenting examples of processes with up to six external legs attached to the loop.


Introduction
The Standard Model is currently being re-discovered at the LHC, and new exclusion limits on Beyond the Standard Model particles-and on the Higgs mass-are being delivered by the experimental collaborations with an impressive speed. Higher order corrections play an important role in obtaining bounds on the Higgs boson and New Physics. In particular, the exclusion limits for the Higgs boson would look very different if we only had leading order tools at hand. Further, it will be very important to have precise theory predictions to constrain model parameters once a signal of New Physics has been established. Therefore it is of major importance to provide tools for next-to-leading order (NLO) predictions which are largely automated, such that signal and background rates for a multitude of processes can be estimated reliably.
We have recently experienced major advances in the activity of constructing packages for fully automated oneloop calculations, see e.g. [38][39][40][41][42][43]. The concepts that lead to these advances have been recently reviewed in [44]. Among the most important developments are the integrand-reduction technique [45,46] and the generalized n-dimensional unitarity [47]. Their main outcome is a numerical reconstruction of a representation of the tensor structure of any one-loop integrand where the multi-particle pole configuration is manifest. As a consequence, decomposing one-loop amplitudes in terms of basic integrals becomes equivalent to reconstructing the polynomial forms of the residues to all multi-particle cuts. Within this algorithm, the integrand of a given scattering amplitude, carrying complete and explicit information on the chosen dimensionalregularisation scheme, is the only input required to accomplish the task of its evaluation. In fact, the integration is substituted by a much simpler operation, namely by polynomial fitting, which requires the sampling of the integrand on the solutions of generalised on-shell conditions.
In this article, we present the program package GOSAM which allows the automated calculation of one-loop amplitudes for multi-particle processes. Amplitudes are expressed in terms of Feynman diagrams, where the integrand is generated analytically using QGRAF [2], FORM [48], spinney [49] and haggies [50]. The individual program tasks are steered via python scripts, while the user only needs to edit an "input card" to specify the details of the process to be calculated, and launch the generation of the source code and its compilation, without having to worry about internal details of the code generation.
The program offers the option to use different reduction techniques: either the unitarity-based integrand reduction as implemented in SAMURAI [40] or traditional tensor reduction as implemented in Golem95C [51,52] interfaced through tensorial reconstruction at the integrand level [53], or a combination of both. It can be used to calculate oneloop corrections within both QCD and electroweak theory. Beyond the Standard Model theories can be interfaced using FeynRules [54] or LanHEP [55]. The Binoth-Les Houchesinterface [56] to programs providing the real radiation contributions is also included.
The advantage of generating analytic expressions for the integrand of each diagram gives the user the flexibility to organize the computation according to his own efficiency preferences. For instance, the computing algorithm can proceed either diagram-by-diagram or by grouping diagrams that share a common set of denominators (suitable for a unitarity-based reduction), and it can deal with the evaluation of the rational terms either on the same footing as the rest of the amplitude, or through an independent routine which evaluates them analytically. These options and the other features of GOSAM will be discussed in detail in the following.
In Sect. 2, after giving an overview on the diagram generation and on processing gauge-group and Lorentz algebra, we discuss the code generation and the reduction strategies.
The installation requirements are given in Sect. 3, while Sect. 4 describes the usage of GOSAM, containing all the set-up options which can be activated by editing the input card. In Sect. 5 we show results for processes of various complexity. The release of GOSAM is accompanied by the generated code for some example processes, listed in Appendix A.

Overview and algorithms
2.1 Overview GOSAM produces, in a fully automated way, all the code required to perform the calculation of one-loop matrix elements. There are three main steps in the process of constructing the code: the generation of all contributing diagrams within a process directory, the generation of the Fortran code, and finally compiling and linking the generated code. These steps are self-contained in the sense that after each step all the files contained in the process directory could be transfered to a different machine where the next step will be carried out.
In the following sections we focus on the algorithms that are employed for the construction of the code to produce and evaluate matrix elements.
The first step (setting up a process directory), which consists in the generation of some general source files and the generation of the diagrams, is described in Sect. 2.2. The second step (generating the fortran code) is carried out by means of advanced algorithms for algebraic manipulation and code optimization which are presented in Sects. 2.3 and 2.4. The third step (compilation and linking) is not specific to our code generation, therefore will not be described here.
The practical procedures to be followed by the user in generating the code will be given in Sect. 4, which can be considered a short version of the user manual.

Generation and organisation of the diagrams
For the diagram generation both at tree level and one-loop level we employ the program QGRAF [2]. This program already offers several ways of excluding unwanted diagrams, for example by requesting a certain number of propagators or vertices of a certain type or by specifying topological properties such as the presence of tadpoles or on-shell propagators. Although QGRAF is a very reliable and fast generator, we extend its possibilities by adding another level of analysing and filtering over diagrams by means of Python. This gives several advantages: first of all, the possibilities offered by QGRAF are not always sufficient to distinguish certain classes of diagrams (see examples in Fig. 1); secondly,  Fig. 1(a) is zero in dimensional regularisation. However, in QGRAF there is no operator to identify this type of diagrams. In Fig. 1(b) the Z boson is emitted from a closed quark line. These diagrams form a separate gauge invariant class and could be treated separately from diagrams where the Z boson comes from an external quark line QGRAF cannot handle the sign for diagrams with Majorana fermions in a reliable way; finally, in order to fully optimize the reduction, we want to classify and group diagrams according to the sets of their propagators.
Within our framework, QGRAF generates three sets of output files: an expression for each diagram to be processed with FORM [48], Python code for drawing all diagrams, and Python code for computing the properties of each diagram. The information about the model for QGRAF is either read from the built-in Standard Model file or is generated from a user defined LanHEP [55] or Universal FeynRules Output (UFO) [54] file.
The Python program automatically performs several operations: -diagrams whose color factor turns out to be zero are dropped automatically; -the fermion flow is determined and used to compute an overall sign for each diagram, which is relevant in the presence of Majorana fermions; -the number of propagators containing the loop momentum, i.e. the loop size of the diagram, the tensor rank and the kinematic invariants of the associated loop integral are computed; -diagrams with an associated vanishing loop integral (see Fig. 1(a)) are detected and flagged for the diagram selection; -all propagators and vertices are classified for the diagram selection; diagrams containing massive quark self-energy insertions or closed massless quark loops are specially flagged.
Any one-loop diagram can be written in the form where the numerator is a polynomial of tensor 1 rank r.
..μ r r q μ 1 · · · q μ r , (2) and the N × N kinematic matrix is defined as All masses can be either real or complex. Important information about the integrals that will appear in the reduction of each one-loop diagram is contained in the tensor rank r of the loop integral and its kinematic matrix S ij . We define a preorder relation on one-loop diagrams, such that D 1 D 2 if their associated matrices S(D 1 ) and S(D 2 ) are related by a finite (not necessarily unique) chain of transformations where each transformation is one of the following: -the identity, -the simultaneous permutation of rows and columns, -the simultaneous deletion of the row and column with the same index, which corresponds to pinching the corresponding propagator in the diagram.
The relation " " can be read as "appears in the reduction of". Our algorithm groups the one-loop diagrams The partitioning procedure provides an important gain in efficiency, because while carrying out the tensor reduction for the diagram max V i , all other diagrams in the same cell V i are reduced with virtually no additional computational cost. The gain in efficiency can be observed when reducing the diagram using the OPP method [45] and its implementations in CutTools [57] and SAMURAI [40], as well as in classical tensor reduction methods as implemented e.g. in Go-lem95C [51,52], PJFRY [58] and LoopTools [3,59]. In order to draw the diagrams, we first compute an ordering of the external legs which allows for a planar embedding of the graph. Such ordering can always be found for a tree or a one-loop graph since non-planar graphs only start to appear in diagrams with two or more loops. After the legs have been assigned to the vertices of a regular polygon, we use our own implementation of the algorithms described in [60] for fixing the coordinates of the remaining vertices; the algorithm has been extended to determine an appealing layout also for graphs containing tadpoles. Starting from these 1 Index contractions in (2) are understood in n-dimensional space. coordinates and using the package Axodraw [61], GOSAM generates a L A T E X file that contains graphical representations of all diagrams.

Color algebra
In the models used by GOSAM, we allow one unbroken gauge group SU(N C ) to be treated implicitly; any additional gauge group, broken or unbroken, needs to be expanded explicitly. Any particle of the model may be charged under the SU(N C ) group in the trivial, (anti-)fundamental or adjoint representation. Other representations are currently not implemented.
For a given process we project each Feynman diagram onto a color basis consisting of strings of generators i p−1 j and Kronecker deltas δ ij but no contractions of adjoint indices and no structure constants f ABC . Considering, for example, the process GOSAM finds the color basis are the color parts of the quark and gluon wave functions respectively. The dimension of this color basis for N g external gluons and N qq quark-antiquark pairs is given by [62]: It should be noted that the color basis constructed in this way is not a basis in the mathematical sense, as one can find linear relations between the vectors |c i once the number of external partons is large enough. Any Feynman diagram can be reduced to the form for the process specific color basis |c 1 , . . . , |c k by applying the following set of relations: The same set of simplifications is used to compute the matrices c i |c j and c i |T I · T J |c j . The former is needed for squaring the matrix element, whereas the latter is used to provide color correlated Born matrix elements which we use for checking the IR poles of the virtual amplitude and also to provide the relevant information for parton showers like POWHEG [33,34,63]. For the above example, GOSAM obtains 2 Similarly, the program computes the matrices c i |T I · T J |c j for all pairs of partons I and J . If M (0) denotes the tree-level matrix element of the process and we have then the square of the tree level amplitude can be written as For the interference term between leading and next-toleading order we use a slightly different philosophy. First of all we note that it is sufficient to focus on a single group V α as defined in Sect. 2.2, + h.c. (12) In order to reduce the complexity at the level of the reduction, we perform the contraction with the tree-level already at the integrand level, where C (1) j is formed by the sum over the corresponding coefficients of all diagrams D ∈ V α .

Lorentz algebra
In this section we discuss the algorithms used by GOSAM to transform the coefficients C (0) i and C (1) i (q), as defined in the previous section, such that the result is suitable for efficient numerical evaluation. One of the major goals is to split 2 In the actual code the results are given in terms of T R and N C only. the n-dimensional algebra (n = 4 − 2ε) into strictly fourdimensional objects and symbols representing the higherdimensional remainder.
In GOSAM we have implemented the 't Hooft-Veltman scheme (HV) and dimensional reduction (DRED). In both schemes all external vectors (momenta and polarisation vectors) are kept in four dimensions. Internal vectors, however, are kept in the n-dimensional vector space. We adopt the conventions used in [49], wherek denotes the fourdimensional projection of an in general n-dimensional vector k. The (n − 4)-dimensional orthogonal projection is denoted ask. For the integration momentum q we introduce in addition the symbol μ 2 = −q 2 , such that We also introduce suitable projectors by splitting the metric tensor g μν =ĝ μν +g μν ,ĝ μνg νρ = 0, In the following, we describe the 't Hooft algebra in detail. For DRED, the only differences are that the numerator algebra is performed in four dimensions for both external and internal vectors (i.e. q ≡q) and that in the very end all appearances of q 2 are replaced byq 2 − μ 2 .
Wave functions and propagators GOSAM contains a library of representations of wave functions and propagators up to spin two. 3 The exact form of the interaction vertices is taken from the model files.
The representation of all wave functions with non-trivial spin is based on massless spinors. Each massive external vector p i is replaced by its light-cone projection l i with respect to a lightlike reference vector k, For spin 1 / 2 particles we use the assignment of wave functions as shown in Table 1; here, we quote the definition of the massive spinors from [49] assuming the splitting of (16): Assignment of quark and lepton wave functions. We label the physical spin states by j 3 = ±1, which is twice the 3-component of the spin. The wave functions assigned in Table (a) are mapped onto the bracket notation used in spinney [49] as defined in Tables (b) and (c) (a) Assignment of initial and final states for quarks and leptons Wave functions for massless fermions In order to preserve the condition that for any loop integral the tensor rank does not exceed the number of loop propagators we fix all gauge boson propagators to be in Feynman gauge. Their wave functions are constructed as [64] ε μ (p, +1) = where p = p in the massless case and p = l according to (16) in the massive case. In the latter case the third polarisation is defined as The wave functions and propagators for spin 3 / 2 and spin 2 particles correspond to those in [65].
Simplifications Once all wave functions and propagators have been substituted by the above definitions and all vertices have been replaced by their corresponding expressions from the model file then all vector-like quantities and all metric tensors are split into their four-dimensional and their orthogonal part. As we use the 't Hooft algebra, γ 5 is defined as a purely four-dimensional object, γ 5 = i μνρσγ μγ νγ ργ σ . By applying the usual anti-commutation relations for Dirac matrices we can always separate the fourdimensional and (n − 4)-dimensional parts of Dirac traces, as we can use the fact that [49,62] tr(1) · tr(γ μ 1 · · ·γ μ lγ μ l+1 · · ·γ μ l+p ) = tr(γ μ 1 · · ·γ μ l ) · tr(γ μ l+1 · · ·γ μ l+p ).
While the (n − 4)-dimensional traces are reduced completely to products of (n − 4)-dimensional metric tensors g μν , the four-dimensional part is treated such that the number of terms in the resulting expression is kept as small as possible. Any spinor line or trace is broken up at any position where a light-like vector appears. Furthermore, Chisholm identities are used to resolve Lorentz contractions between both Dirac traces and open spinor lines. If any traces remain we use the built-in trace algorithm of FORM [48].
In the final result we can always avoid the explicit appearance of Levi-Civitá tensors, noticing that any such tensor is contracted with at least one light-like vector 4kμ , and we can replacê Hence, the kinematic part of the numerator, at the end of our simplification algorithm, is expressed entirely in terms of: -spinor products of the form k i k j , [k i k j ] or [k i |γ μ |k j · q μ , -dot productsk i ·k j ork i ·q, -constants of the Lagrangian such as masses, widths and coupling constants, -the symbols μ 2 =q 2 − q 2 and ε = (n − 4)/2.
Treatment of R 2 rational terms In our representation for the numerator of one-loop diagrams, terms containing the symbols μ 2 or ε can lead to a so-called R 2 term [66], which contributes to the rational part of the amplitude. In general, there are two ways of splitting the numerator function: or, alternatively, It should be noted that in (23a) the terms N 1 and N 2 do not arise in DRED, where only terms containing μ 2 contribute to R 2 . Instead of relying on the construction of R 2 from specialized Feynman rules [67][68][69][70], we generate the R 2 part along with all other contributions without the need to separate the different parts. For efficiency reasons, however, we provide an implicit and an explicit construction of the R 2 terms. The implicit construction uses the splitting of (23a) and treats all three numerator functions N i on equal grounds. Each of the three terms is reduced separately in a numerical reduction and the Laurent series of the three results are added up taking into account the powers of ε.
The explicit construction of R 2 is based on the assumption that each term inÑ in (23b) contains at least one power of μ 2 or ε. The expressions for those integrals are relatively simple and known explicitly. Hence, the part of the amplitude which originates fromÑ is computed analytically whereas the purely four-dimensional partN is passed to the numerical reduction.

Abbreviation system
To prepare the numerator functions of the one-loop diagrams for their numerical evaluation, we separate the symbol μ 2 and dot products involving the momentumq from all other factors. All subexpressions which do not depend on either q or μ 2 are substituted by abbreviation symbols, which are evaluated only once per phase space point. Each of the two parts is then processed using haggies [50], which generates optimized Fortran code for their numerical evaluation. For each diagram we generate an interface to SAMU-RAI [40], Golem95C [52] and/or PJFRY [58]. The two latter codes are interfaced using tensorial reconstruction at the integrand level [53].

Reduction strategies
In the implementation of GOSAM, great emphasis has been put on maintaining flexibility with respect to the reduction algorithm that the user decides to use. On the one hand, this is important because the best choice of the reduction method in terms of speed and numerical stability can strongly depend on the specific process. On the other hand, we tried to keep the code flexible to allow further extensions to new reduction libraries, such that GOSAM can be used as a laboratory for interfacing future methods with a realistic environment.
Our standard choice for the reduction is SAMURAI, which provides a very fast and stable reduction in a large part of the phase space. Furthermore, SAMURAI reports to the client code if the quality of the reconstruction of the numerator suffices the numerical requirements (for details we refer to [40]). In GOSAM we use this information to trigger an alternative reduction with either Golem95C [52] or PJFRY [58] whenever these reconstruction tests fail, as shown in Fig. 2. The reduction algorithms implemented in these libraries extend to phase space regions of small Gram determinants and therefore cover most cases in which onshell methods cannot operate sufficiently well. This combination of on-shell techniques and traditional tensor reduction is achieved using tensorial reconstruction at the integrand level [53], which also provides the possibility of running on-shell methods with a reconstructed numerator. In addition to solving the problem of numerical instabilities, in some cases this option can reduce the computational cost of the reduction. Since the reconstructed numerator is typically of a form where kinematics and loop momentum dependence are already separated, the use of a reconstructed numerator tends to be faster than the original procedure, in particular in cases with a large number of legs and low rank.
The flowchart in Fig. 2 summarizes all possible reduction strategies which are currently implemented. The strategy in use is selected by assigning the variable reduc-tion_interoperation in the generated Fortran code. The availability of the branches is determined during code generation by activating (at least one of) the extensions (samurai, golem95, pjfry) in the input card. SAMURAI with rescue option (Golem95C); 3: SAMURAI with numerator from tensorial reconstruction; 4: same as 3 but with rescue option (Golem95C). 11, 12 and 14 are the same as 1, 2, 3 (respectively) with the difference that PJFRY is used instead of Golem95C. Switching between active branches is possible at run time. In detail, the possible choices for the variable reduc-tion_interoperation are the following: 0 the numerators of the one-loop diagrams are reduced by SAMURAI, no rescue system is used in case the reconstruction test fails; 1 the tensor coefficients of the numerators are reconstructed using the tensorial reconstruction at the integrand level, the numerator is expressed in terms of tensor integral form factors which are evaluated using Go-lem95C; 2 the numerators are reduced by SAMURAI; whenever the reconstruction test fails, numerators are reduced using the option 1 as a backup method; 3 tensorial reconstruction is used to compute the tensor coefficients; SAMURAI is employed for the reduction of the reconstructed numerator, no rescue system is used; 4 as in option 3, SAMURAI is used to reduce the reconstructed numerator, Golem95C is used as backup option; 11 same as 1 but PJFRY is used instead of Golem95C; 12 same as 2 but PJFRY is used instead of Golem95C; 14 same as 4 but PJFRY is used instead of Golem95C.
It is difficult to make a statement about the "optimal" reduction method because this depends on the process under consideration. For multi-leg processes, e.g. bbbb production, we found that SAMURAI is clearly superior to tensor reduction in what concerns timings and size of the code. Concerning points which need a special treatment, we did not make extensive studies using traditional tensor reduction only, but one can certainly say that the combination of SA-MURAI and tensorial reconstruction seems to be optimal in what concerns the avoidance of numerical instabilities due to inverse Gram determinants.

Conventions of the amplitudes
In this section we briefly discuss the conventions chosen for the results returned by GOSAM. Depending on the actual setup for a given process, in particular if an imported model file is used, conventions may be slightly different. Here we restrict the discussion to the case where the user wants to compute QCD corrections to a process and in the setup files he has put g s = 1. In this case, the tree-level matrix element squared can be written as The fully renormalised matrix element at one-loop level, i.e. the interference term between tree-level and one-loop, can be written as A call to the subroutine samplitude returns an array consisting of the four numbers (a 0 , c 0 , c −1 , c −2 ) in this order. The average over initial state colours and helicities is included in the default setup. In cases where the process is loop induced, i.e. the tree level amplitude is absent, the program returns the values for has been pulled out.
After all UV-renormalisation contributions have been taken into account correctly, only IR-singularities remain, which can be computed using the routine ir_subtractions. This routine returns a vector of length two, containing the coefficients of the single and the double pole, which should be equal to (c −1 , c −2 ) and therefore can be used as a check of the result.
Ultraviolet renormalisation in QCD For UV-renormalisation we use the MS scheme for the gluon and all massless quarks, whereas a subtraction at zero momentum is chosen for massive quarks [71]. Currently, counterterms are only provided for QCD corrections. In the case of electroweak corrections only unrenormalised results can be produced automatically.
For computations involving loop propagators for massive fermions, we introduced the automatic generation of a mass counter term needed for the on-shell renormalisation of the massive particle. Here, we exploit the fact that such a counter term is strictly related to the massive fermion self energy bubble diagrams (see Fig. 3). As described in Sect. 2.2, the program GOSAM analyzes all generated diagrams. In that step also self-energy insertions of massive The symbol 1 HV is one in the 't Hooft Veltman scheme and zero in DRED.
Performing the integral, contracting the expression with the QCD vertices at both sides and multiplying the missing factor of (2π) −1 we retrieve the expression for the mass counter-term, Furthermore, the renormalisation of α s leads to a term of the form with β 0 = (11C A − 4T R N f )/6, N f being the number of light quark flavours, N f,h the number of heavy flavours, and b is the power of the coupling in the Born amplitude as defined in (24). The last term of (28) provides the finite renormalisation needed to compensate the scheme dependence of α s , A further contribution consists of the wave-function renormalisation of massive external quark lines. If we denote the set of external massive quark lines by Finally, also the wave function of the gluon receives a contribution from the presence of heavy quarks in closed fermion loops. If N g is the number of external gluon lines, this contribution can be written as At the level of the generated Fortran code the presence of these contributions can be controlled by a set of variables defined in the module config.f90. The variable renormalisation can be set to 0, 1, or 2. If renormalisa-tion=0, none of the counterterms are present. If renor-malisation=2 only |M| 2 ct,δm Q is included, which is the counterterm stemming from all terms of the type of (27) contributing to the amplitude. In the case where renormalisation=1 a more finegrained control over the counterterms is possible.
renorm_logs: if set to false, in all counterterms the generation of logarithms is disabled, i.e. factors of the form (•) ε in (27) to (31) are replaced by one. renorm_beta: if set to false, the counterterm |M| 2 ct,α s is set to zero. renorm_mqwf: if set to false, the counterterm |M| 2 wf,Q is set to zero. renorm_mqse: if set to false, the counterterm |M| 2 ct,δm Q is set to zero. renorm_decoupling if set to false, the counterterm |M| 2 wf,g is set to zero. The default settings for renormalisation=1 are true for all the renorm options listed above.
Finite renormalisation of γ 5 in QCD In the 't Hooft Veltman scheme, a finite renormalisation term for γ 5 is required beyond tree level. The relevant terms are generated only if fr5 is added in the input card to the list of extensions before code generation. Currently, the automatic generation of this finite contribution is not performed if model files different from the built-in model files are used. In agreement with [72] and [73] we replace the axial component at each vertex, with Once it is generated, this contribution can be switched on and off at run-time through the variable renorm_gamma5, which is defined in the module config.f90.
Conversion between the schemes In GOSAM we have implemented two different schemes, the 't Hooft Veltman scheme and dimensional reduction. By default, the former is used, while the latter can be activated by adding the extension dred. If a QCD computation has been done in dimensional reduction the result can be converted back to the 't Hooft Veltman scheme by adding a contribution for each external massless parton, . This conversion can be switched on by setting convert_to_cdr to true in the module config.f90. At one-loop level, the 't Hooft Veltman scheme and conventional dimensional regularisation (CDR) are equivalent in the sense thatγ t HV I = 0 for all partons.

Requirements
The program GOSAM is designed to run in any modern Linux/Unix environment; we expect that Python (≥ 2.6), Java (≥1.5) and Make are installed on the system. Furthermore, a Fortran 95 compiler is required in order to compile the generated code. Some Fortran 2003 features are used if one wants to make use of the Les Houches interface [56]. We have tried all examples using gfortran versions 4.1 and 4.5.
On top of a standard Linux environment, the programs FORM [48], version ≥ 3.3 (newer than Aug. 11, 2010) and QGRAF [2] need to be installed on the system. Whereas spinney [49] and haggies [50] are part of GOSAM and are not required to be installed separately, at least one of the libraries SAMURAI [40] and Golem95C [52] needs to be present at compile time of the generated code. Optionally, PJFRY [58] can be used on top of Golem95C.

Download and installation
QGRAF The program can be downloaded as Fortran source code from http://cfif.ist.utl.pt/~paulo/qgraf.html.
After unpacking the tar-ball, a single Fortran 77 file needs to be compiled. The build process and installation of GOSAM is controlled by Python Distutils, while the build process for the libraries SAMURAI and Golem95C is controlled by Autotools. Therefore the installation proceeds in two steps: 1. For all components which use Autotools, the following sequence of commands installs them under the user defined directory MYPATH.
./configure --prefix=MYPATH make FC=gfortran F77=gfortran make install # or sudo make install If the configure script is not present, the user needs to run sh ./autogen.sh first. 2. For GOSAM which is built using Distutils, the user needs to run If MYPATH is different from the system default (e.g. /usr/bin), the environment variables PATH, LD_LIBRA-RY_PATH and PYTHONPATH might have to be set accordingly. For more details we direct the user to the GO-SAM reference manual and to the documentation of the beforementioned programs.

Setting up a simple process
GOSAM is a very flexible program and comes with a wide range of configuration options. Not all of these options are relevant for simple processes and often the user can leave most of the settings at their default values. In order to generate the code for a process, one needs to prepare an input file, which will be called process card in the following, which contains -process specific information, such as a list of initial and final state particles, their helicities (optional) and the order of the coupling constants; -scheme specific information and approximations, such as the regularisation and renormalisation schemes, the underlying model, masses and widths which are set to zero, the selection of subsets of diagrams; the latter might be process dependent; -system specific information, such as paths to programs and libraries or compiler options; -optional information for optimisations which control the code generation.
In the following we explain how to set up the required files for the process qq → gZ 0 → g e − e + . The example computes the QCD corrections for the uū initial state, where m e = 0 and N f = 5 massless quarks are assumed. For our example, we follow an approach where we keep the different types of information in separate files-process.rc, scheme.rc and system.rc-and use GOSAM to produce a process card for this process based on these files. This is not required-one could also produce and edit the process card directly-it is however more convenient to store system specific information into a separate, re-usable file, and it makes the code generation more transparent.

Process specific information
The following listing contains the information which is specific to the process. The syntax of process cards requires that no blank character is left between the equals sign and the property name. Commentary can be added to any line, marked by the '#' character. Line continuation is achieved using a backslash at the end of a line. 5 Listing 1 File 'process.rc' 1 process_path=qqgz 2 in=u,u3 out=g,e-,e+ 4 helicities=+-+-+,+---+,-++-+,-+--+ 5 order=QCD,1,3 5 The line numbers are just for reference and should not be included in the actual files.
The first line defines the (relative) path to the directory where the process files will be generated. GOSAM expects that this directory has already been created. Lines 2 and 3 define the initial and final state of the process in terms of field names, which are defined in the model file. Besides the field names one can also use PDG codes [76,77] instead. Hence, the following lines would be equivalent to lines 2 and 3 in Listing 1: 2 in=2,-2 3 out=21,11,-11 Line 4 describes the helicity amplitudes which should be generated. If no helicities are specified, the program defaults to the generation of all possible helicity configurations, some of which may turn out to be zero. The different helicity amplitudes are separated by commas; within one helicity amplitude there is one character (usually '+', '-' and '0') per external particle from the left to the right. In the above example for the reaction we have the following assignments: With the above value for helicities we generate all non-vanishing helicities for the partons but keep the lepton helicities fixed. In more complicated examples this way of listing all helicities explicitly can be very tedious. Therefore, we introduced the option to generate sets of helicities using square brackets. For example, if the gluon helicity is replaced by [+-], the bracket is expanded automatically to take the values +,-. Finally, the order (power) of the coupling constants has to be specified. Line 5 contains a keyword for the type of coupling (QCD or QED), the order of this coupling constant in the unsquared tree level amplitude (in our example: 1) and the order of the coupling constant in the unsquared one-loop amplitude (in our example: 3). One can also use GOSAM to generate the tree level only, by giving only the power of the tree level amplitude: Conversely, GOSAM will generate the virtual amplitude squared for processes where no tree level is present if the tree level order is replaced by the keyword NONE. The alternative solution is the application of a Python filter using the options filter.lo for tree level and filter.nlo for one-loop diagrams. The current example requires the two lines In Listing 2, line 1 selects dimensional reduction as a regularisation scheme. If dred is not specified in the list of extensions, GOSAM works in the 't Hooft Veltman scheme by default. Line 2 removes all on-shell bubbles on external legs. This is, on the one hand, required to be consistent with the renormalisation scheme. On the other hand, those diagrams would lead to zero denominators at the algebraic level. In line 3 all light quark masses, the mass of the electron and the width of the top quark are set to zero. Further, as a convention rather than a scheme, the strong coupling g s is set to one in line 4, which means that g s will not occur in the algebraic expressions, assuming that the user will multiply the final result by his desired value for the strong coupling. If the option one=gs is not used, the default value contained in the file common/model.f90 will be used. This default value of course can be changed by the user.
System specific information In order to adapt the code generation to the system environment, GOSAM needs to find a way of determining all relevant paths and options for the programs and libraries used during generation, compilation and linking of the code. Those settings are fixed in the file system.rc in our example. Generating the code After having prepared the input files correctly we need to collect the information distributed over the three files process.rc, scheme.rc and system.rc in one input file, which we will call gosam.in here. The corresponding command is: gosam.py --template gosam.in \ --merge process.rc \ --merge scheme.rc --merge system.rc The generated file can be processed with gosam.py directly but requires the process directory to be present.
mkdir qqgz gosam.py gosam.in cd qqgz All further steps are controlled by the generated make files; in order to generate and compile all files relevant for the matrix element one needs to invoke make compile The generated code can be tested with the program matrix/test.f90. The following sequence of commands will compile and run the example program. The printed numbers are, in this order, a 0 , c 0 /a 0 , c −1 /a 0 , c −2 /a 0 and the pole parts calculated from the infrared insertion operator [78,79].
One can generate a visual representation of all generated diagrams using the command make doc which generates the file doc/process.ps using a Python implementation of the algorithm described in [60] and the L A T E X package AXODRAW [61].

Further options
GOSAM provides a range of options which influence the code generation, the compilation and the numerical evaluation of the amplitude. Giving an exhaustive list of all options would be far beyond the scope of this article and the interested user is referred to the reference manual. Nonetheless, we would like to point out some of GOSAM's capabilities by presenting the corresponding options. 7 The actual numbers depend on the random number generator of the system because the phase space point is generated randomly; however, the pole parts should agree between the matrix element and the infrared insertion operator given that the matrix element is fully renormalised.
Generating the R 2 term When setting up a process the user can specify if and how the R 2 term of the amplitude should be generated by setting the variable r2 in the setup file.

r2=explicit
Possible options for r2 are implicit, which is the default, explicit, off and only. The keyword implicit means that the R 2 term is generated along with the four-dimensional numerator as a function in terms ofq, μ 2 and ε and is reduced at runtime by sampling different values for μ 2 . This is the slowest but also the most general option. Using the keyword explicit carries out the reduction of terms containing μ 2 or ε during code generation (see Appendix B). The keyword off puts the R 2 term to zero which is useful if the user wants to provide his own calculation for these terms. Conversely, using r2=only discards everything but the R 2 term (reducing it as in the case explicit) and puts GOSAM in the position of providing R 2 terms for external codes which work entirely in four dimensions.
Diagram selection GOSAM offers a two-fold way of selecting and discarding diagrams. One can either influence the way QGRAF generates diagrams or apply filters to the diagrams after they have been generated by QGRAF or combine the two methods. Let us assume that in the above example we want to remove the third generation of quarks completely. Hence, all closed quark loops would be massless and therefore the second generation is just an exact copy of the first one. We can therefore restrict the generation of closed quark loops to up and down quarks. GOSAM has a filter precisely for this purpose, which takes the field names of the flavours to be generated as arguments.

filter.nlo=NFGEN(U,D)
This filter can be combined with the already existing filter selecting only diagrams containing a Z-propagator using the AND function:

filter.nlo=AND( NFGEN(U,D), \ IPROP([Z]) == 1 )
A further feature of the code generated by GOSAM is the possibility of selecting diagrams at runtime. For example, we would like to distinguish at runtime three different gauge invariant sets of diagrams at one-loop level: 1. diagrams with a closed quark loop where the Z is attached to the loop; 2. diagrams with a closed quark loop where the Z is emitted from the external quark line; 3. diagrams without a closed quark loop.
In order to provide the code for a diagram selection at runtime one simply replaces the above filter by a list of filters as follows Additional extensions Some of GOSAM's functionality is available through the extensions variable. On top of the already presented options for selecting a regularisation scheme (by adding the option dred) or for activating interfaces to several different reduction libraries (samurai, golem95, pjfry) the user can also add the following options: fr5 adds and activates the relevant code for the computation of the finite renormalisation of γ 5 required in the 't Hooft Veltman scheme as described in (32). powhegbox generates routines for the computation of the color and spin correlated Born matrix elements as required by POWHEG [34]. autotools uses make files which use Autoconf and Automake for compilation of the matrix element. gaugecheck replaces the polarisation vectors of external vector fields by where the variable z i defaults to zero and is accessible in the Fortran code through the symbol gaugeiz.

Interfacing the code
The matrix element code generated by GOSAM provides several routines to transparently access partial or full results of the amplitude calculation. Here, we only present a minimal set of routines which can be used to obtain the set of coefficients [a 0 , c 0 , c −1 , c −2 ] for a given scale and a given set of external momenta. The routines, which can be accessed through the modules matrix 8 are defined as follows: initgolem This subroutine must be called once before the first matrix element evaluation. It initializes all dependent model parameters and calls the initialisation routines of the reduction libraries. The first mandatory arguments of this routine are the external momenta vecs, where vecs(i,:) contains the momentum of the i-th particle as a vector [E i , p x i , p y i , p z i ], and we use in-out kinematics, i.e. p 1 + p 2 → p 3 + · · · + p N . Maximal numerical stability is achieved if the beam axis is chosen along the z-axis. The second argument, scale2 = μ 2 R , is the square of the renormalisation scale. As a third argument the routine expects a vector which accepts the result in the format [a 0 , c 0 , c −1 , c −2 ] with the coefficients being defined in (24) and (25). The optional argument ok may be used in order to report the outcome of the reconstruction tests in samurai if no rescue method has been chosen (see Sect. The optional argument exit_libs should only be set if multiple calls to this routine (e.g. for different matrix elements) are necessary and the dependent libraries should be terminated only once.
A small program which computes the amplitude for a set of phase space points is automatically generated with the amplitude code in the file test.f90 in the subdirectory matrix. The script config.sh in the process directory returns suitable compilation and linking options for the generated matrix element code.

Using the BLHA interface
The so-called Binoth Les Houches Accord (BLHA) [56] defines an interface for a standardized communication between one-loop programs (OLP) and Monte Carlo (MC) tools. The communication between the two sides is split into two main phases: an initialisation phase and a runtime phase. During initialisation the two programs establish an agreement by exchanging a set of files and typically initiate the code generation. The OLP runtime code is then linked to the MC program and, during the runtime phase, called through a welldefined set of routines providing NLO results for the phase space points generated by the MC. According to this standard, it is the responsibility of the MC program to provide results for the Born matrix element, for the real emission and for a suitable set of infrared subtraction terms. A schematic overview on this procedure is shown in Fig. 4.
GOSAM can act as an OLP in the framework of the BLHA. In the simplest case, the MC writes an order file-in this example it is called olp_order.lh-and invokes the script gosam.py as follows: gosam.py --olp olp_order.lh Further, GOSAM specific options can be passed either in a file or directly at the command line. One can, for example, use autotools for the compilation by modifying the above line as follows.

gosam.py --olp olp_order.lh \ extensions=autotools
The contract file is given the extension .olc by default and would be olp_order.olc in this example. Alternatively, the name can be altered using the -o option.
If successful, the invocation of gosam.py generates a set of files which can be compiled as before with a generated make file. The BLHA routines are defined in the Fortran module olp_module but can also be accessed from C programs. 9 The routines OLP_Start and OLP_EvalSubProcess are defined exactly as in the BLHA proposal [56]. For convenience, we extended the interface by the functions OLP_Finalize(), which terminates all reduction libraries, and OLP_Option(char*, int*), which can be used to pass non-standard options at runtime. For example, a valid call in C to adjust the Higgs mass would be int ierr; OLP_Options("mH=146.78", &ierr); A value of one in ierr indicates that the setting was successful. A value of zero indicates an error.

Using external model files
With a few modifications in the process description files, GOSAM can immediately make use of model files generated by either FeynRules [80] in the UFO format [54] or by LanHEP [55]. In both cases, the following limitations and differences with respect to the default model files, sm and smdiag, apply: are not guaranteed to work with models other than the built-in models.
Importing models in the UFO format A model description in the UFO format consists of a Python package stored in a directory. In order to import the model into GOSAM one needs to set the model variable in the input card to specify the keyword FeynRules in front of the directory name, where we assume that the model description is in the directory $HOME/models/MSSM_UFO.

model=FeynRules,$HOME/models/MSSM_UFO
Importing models in the LanHEP format LanHEP model descriptions consist of a set of plain text files in the same directory with a common numbering (such as func4.mdl, lgrng4.mdl, prtcls4.mdl, vars4.mdl). A Lan-HEP model can be loaded by specifying the path and the common number in the model variable. Assuming the files are situated in the directory $HOME/models/ MSSM_LHEP one would set the variable as follows.

model=$HOME/models/MSSM_LHEP,4
Details about the allowed names for the table columns are described in the GOSAM reference manual. Precompiled MSSM_UFO and MSSM_LHEP files can also be found in the subdirectory examples/model.

Sample calculations and benchmarks
The codes produced by GOSAM have been tested on several processes. In this section we describe some examples of applications. Additional results, whose corresponding code is also included in the official distribution of the program, will be reported in Appendix B.

pp → W − + j with SHERPA
In Sect. 4.3 the BLHA interface of GOSAM was presented. This interface allows one to link the program to a Monte Carlo event generator, which is, in general, responsible for supplying the missing ingredients for a complete NLO calculation of a physical cross section. Among the different general purpose Monte Carlo event generators, SHERPA [81] is one of those which offers these tools: computing the LO cross section, the real corrections with both the subtraction terms and the corresponding integrated counterparts [82][83][84]. Furthermore, SHERPA offers the possibility to match a NLO calculation with a parton shower [85,86]. Using the BLHA interface, we linked GOSAM with SHERPA to compute the physical cross section for W − +1 jet at NLO. The first steps to perform this linking is to write a SHER-PA input card for the desired process. Instructions and many examples on how to write this can be found in the on-line manual [87]. Running the code for the first time will produce an order file OLE_order.lh which contains all the necessary information for GOSAM, to produce the desired code for the loop part of the process. This includes a list of all partonic subprocesses needed. In parallel to the production of the needed SHERPA libraries with the provided script, one can at this point run the gosam.py command with the flag -olp and the correct path to the order file as explained in Sect. 4.3. Further options may be specified. Among them it is useful to have a second, GOSAM-specific, input card with all the important GOSAM options. Since, at the end, SHER-PA needs to be linked to a dynamic library, it is convenient to run GOSAM with the autotools extension, which allows the direct creation of both static and dynamic libraries, together with the test routine test. The gosam.py script creates all the files needed for interfacing GOSAM with the Monte Carlo event generator together with the code for the one-loop computation of all needed subprocesses, and a makefile to run them. The different parton-level subprocesses are contained in different subdirectories. At this point the user simply has to run the makefile to generate and compile the code. Once the one-loop part of the code is ready, the produced shared library must be added to the list of needed libraries in the SHERPA input card as follows.

SHERPA_LDADD = LHOLE golem_olp;
With this operation the generation of the code is completed. The evaluation of the process and the physical analysis can then be performed at the user's discretion following the advice given in the SHERPA on-line documentation [87].
We tested the BLHA interface by computing W − + 1 jet and producing distributions for several typical observables.
In Figs. 5(a) and 5(b) the inclusive transverse momentum and rapidity of the jets is shown. These distributions were compared with similar ones produced using the program MCFM [27,28], and perfect agreement was found.

pp → W ± + j , EW corrections
As a first example of an electroweak calculation, we computed the virtual one-loop corrections to ud → Wg. A complete analytical calculation for this process was presented in Ref. [ For the kinematic point given in Table 2 (67), (70) of Ref. [88] and with Ref. [89] for the infrared divergences that remain after renormalisation.

γ γ → γ γ
The process γ γ → γ γ in the Standard Model first arises at the one-loop order, and proceeds through a closed loop of fermions and W bosons. Of the 16 helicity amplitudes contributing to it, only three are independent and their analytic expressions can be found in [90]. The pure QED contribution, involving a fermion loop, is contained in samurai-1.0 [40] and will not be repeated here. Instead, we show the results of the W -loop contribution to the independent helicity amplitudes, as an example of EW corrections that can be handled with GOSAM.
With the above parameters and the kinematics of Table 3 we obtain the following results.  [54] which facilitates the import of Feynman rules generated by FeynRules [80] to programs generating one-loop amplitudes. To import such files within the GoSam setup, all the user has to do is to give the path to the corresponding model file in the input card.
For this example, we combined the one-loop amplitude with the real radiation corrections to obtain results for differential cross sections. A calculation of neutralino pair production for the LHC presenting total cross sections at NLO is given in [91].
For the infrared subtraction terms the program Mad-Dipole [92,93] is used, the real emission part is calculated using MadGraph/MadEvent [94]. The virtual matrix element is renormalized in the MS scheme, while massive particles are treated in the on-shell scheme. The renormalisation terms specific to the massive MSSM particles have been added manually.
In Fig. 6 we show the differential cross section for the m χ 0 1 χ 0 1 invariant mass, where we employed a jet veto to suppress large contributions from the channel qg → χ 0 1 χ 0 1 q which opens up at order α 2 α s , but for large p jet T belongs to the distinct process of neutralino pair plus one hard jet production at leading order. We used N f = 5 massless quark flavours and the MSTW08 [95] parton distribution functions. For the SUSY parameters we use the modified benchmarks point SPS1amod suggested in [96], and we use √ s = 7 TeV.
For reference, we also give the result for the unrenormalised amplitude at one specific phase space point for uū → χ 0 1 χ 0 1 in the DRED scheme, using the following parameters and momenta: All widths have been set to zero; for further settings we refer to the model parameter files contained in the subdirectory examples/model/MSSM_UFO. We have checked that the pole terms of the renormalised amplitude cancel with the infrared poles from MadDipole. For the phase  Table 4 we obtain the following numbers.
As an example of a QED calculation, we compared the virtual QED corrections for the process e + e − → e + e − γ with the results provided in [97]. The results compared in the table are the bare unrenormalised amplitudes in the 't Hooft Veltman scheme. No counterterms or subtraction terms have been added to the result.
Using the parameters given above and the kinematics of Table 5 we obtain the following results.

pp → ttH
This process has been compared with the results given in [39]. The partonic subprocesses uū → ttH and gg → ttH where computed both in the 't Hooft Veltman scheme and in dimensional reduction and the fully renormalised results were successfully compared as an internal consistency check. Apart from wave function renormalisation and mass counterterms, Yukawa coupling renormalisation is also needed here. Yukawa coupling counterterms are in this case equal to the wave function counterterms. The Yukawa top mass is set equal to its pole mass.
The kinematics used to obtain the results below is given in Table 6. The results are given in the 't Hooft Veltman scheme, and are fully renormalised. On an Intel Core i7 950 at 3 GHz the evaluation of a single phase space point took 44 ms in the uū channel and 223 ms in the gg channel. The code was compiled with gfortran without optimisations.

gg → ttZ
This amplitude, fully renormalised, has been compared with the results given in [37].
The kinematics used to obtain the results below is given in Table 7.  The evaluation of a single phase space point took 1433 ms on a 2 GHz processor. The code was compiled with gfortran -O2.

pp → bbbb + X
A detailed discussion of this process can be found in [98,99]. In this section we focus on the parts that are relevant in the context of the virtual corrections. In particular we compared our result to the one given in [38], which is the fully renormalised amplitude including the mass counterterms for the top-quark contribution.
The results below are obtained for the phase space point of On an Intel Xeon E7340 the running time for the calculation of a single phase space point was 19.6 s for the gluon initiated channel and 440 ms for the quark channel.

pp → ttbb + X
This process has been compared with the results given in [38]. We have set up the process both in the 't Hooft Veltman scheme and in dimensional reduction and successfully compared the fully renormalised results as an internal consistency check. The results below are given in the 't Hooft Veltman scheme, and only the counterterms for |M| 2 ct,δm t are included.
Using the above parameters and the phase space point of Table 9 we obtain the following results. On an Intel Core i7 950 at 3 GHz the evaluation of a single phase space point took 393 ms in the uū channel and 12.3 s in the gg channel. The code was compiled with gfortran without optimisations.

pp → W + W − bb
The subprocesses uū → W + W − bb and gg → W + W − bb have been calculated both in [38] and [39]. Accordingly, the results below are given in the 't Hooft Veltman scheme, where only the counterterms for |M| 2 ct,δm t are included.
The amplitude ud → W + ggg is an important channel in the calculation of the process pp → W + + 3 jets. The QCD corrections to this process have been presented in Refs. [6][7][8][9].
The subprocess with one quark pair and three gluons consists of more than 1500 Feynman diagrams. We have computed the amplitude including the leptonic decay of the Wboson and compared our result to [38].  On an Intel Core 2 i5 Laptop at 2.0 GHz the evaluation of a single phase space point took about 2.5 s for ud → e + ν e ggg and about 7.5 s for on-shell W's without decay. The code was compiled with gfortran -02.

ud → W + (→ ν e e + )bb (massive b-quark)
The process ud → W + bb, with an on-shell W -boson, has been studied in [100], while the effects of the W -decay have been recently accounted for in [101], and implemented within MCFM. We consider the latter process, and compare the renormalised amplitude evaluated by MCFM. The b-quark is treated as massive in all diagrams except in the vacuum-polarisation like contributions. The evaluation of a single phase space point took 9.12 ms on a 2 GHz processor. The code was compiled with gfortran -O2.

Conclusions
We have presented the program package GOSAM which produces, in a fully automated way, the code required to perform the evaluation of one-loop matrix elements for multiparticle processes. The program is publicly available at http://projects.hepforge.org/gosam/ and can be used to calculate one-loop amplitudes within QCD, electroweak theory, or other models which can be imported via an interface to LanHEP and UFO, also included in the release. Monte Carlo programs for the real radiation can be easily linked through the BLHA interface. GOSAM is extremely flexible, allowing for both unitarity-based reduction at integrand level and traditional tensor reduction, or even for a combination of the two approaches when required. The amplitudes are generated in The calculation of the rational terms is very modular and can proceed either along with the same numerical reduction as the rest of the amplitude, or independently, before any reduction, by using analytic information on the integrals which can potentially give rise to a rational part. In the current version of the code, UV-renormalisation counterterms are provided for QCD corrections only. Further improvements concerning the full automatisation of electroweak corrections are planned.
Different systems to detect and rescue numerical instabilities are implemented, and the user can switch between them without having to re-generate the source code. Due to a careful organisation of the calculation both at the code generation stage and at the reduction stage, the runtimes for multiparticle amplitudes are very satisfactory. Moreover, the GO-SAM generator can also produce codes for processes that include intermediate states with complex masses.
Within the context of the automated matching of Monte Carlo programs to NLO virtual amplitudes, GOSAM can be used as a module to produce differential cross sections for multi-particle processes which can be compared directly to experiment. Therefore we believe that GOSAM can contribute to the goal of using NLO tools as a standard framework for the LHC data analysis at the TeV scale.
Open Access This article is distributed under the terms of the Creative Commons Attribution License which permits any use, distribution, and reproduction in any medium, provided the original author(s) and the source are credited.

Appendix A: Examples included in the release
In the following we give results for the processes listed in the examples directory. Unless stated otherwise, we assume that the coupling constants (e and g s in the standard model) have been set to one in the input card. The conventions for the returned numbers (a 0 , c 0 , c −1 , c −2 ) are as stated in Sect. 2.5. Dimensionful parameters are understood to be in powers of GeV.
As an illustration of the potential of GOSAM, we display in Table 13 the timings required by a wide list of benchmark processes. The first value provided in the table is the time required for the code generation (Generation, given in seconds): we remind the reader that this operation only needs to be performed once per process. The second value is the timing for the full calculation of the amplitude at one phasespace point (Evaluation, in milliseconds). Results are obtained with an Intel(R) Core(TM) i7 CPU 950 @ 3.07 GHz. Table 13 Time required for code generation and calculation of one phase-space point. The results were obtained with an Intel(R) Core(TM) i7 CPU 950 @ 3.07 GHz. The time for the evaluation of a phase space point is taken as the average of the time obtained from the evaluation of 100 random points generated using RAMBO [102], where the code was compiled using gfortran without any optimisation options. The generation of the R 2 term was set to explicit The examples e + e − → tt have an explicit dependence on the Golem95C library and will therefore fail if the extension golem95 is not added.
A.2 e + e − → uu The following parameters and momenta have been used to produce the numerical result: This example has been produced twice: once with the default model file and once with a model file imported from LanHEP [55]. Thus it also can serve as an example of how to import model parameters from LanHEP. The result is given in dimensional reduction, and no renormalisation terms are included. parameters The following results are obtained with the above parameters and the kinematic point of Table 14. This example has been produced twice: once in the 't Hooft Veltman (HV) scheme and once with dimensional reduction (DRED). Only the result in the HV scheme will be Evaluating the amplitude for above parameters and the phase space point given in Table 16   With the above parameters and the kinematics given in Table 19 we obtain the following results. For this process the mass of the b-quark is set to zero. However, in order to have a coupling between the b-quark and the Higgs boson, the following Yukawa coupling is implemented in the model file: With the above parameters and the kinematics given in Table 20 we obtain the following results. The decay width Γ H →γ γ of this loop induced process is known analytically at lowest order. For comparison we used the equations including the top loop and the bosonic contribution given in [108,109]. The decay width can be expressed as A.12 ud → e − ν e g We list the renormalised amplitude in the HV scheme. With the above parameters and the kinematics given in Table 23 we obtain the following results.