This success of gauge theories or standard models of elementary particles, implies that we have definite Lagrangians and thus we can, in principle, predict any process based on these Lagrangians in perturbation theory. When one wants to perform calculation in this way, however, one meets a technical difficulty due to the complexity of the interaction Lagrangian. This is particular to non-abelian gauge theory in which we have three- and four-point self-couplings of gauge bosons as well as interactions of unphysical particles such as Goldstone bosons or ghost particles in general covariant gauge fixing. Hence even in the lowest order of perturbation, that is, in tree level, one finds a number of diagrams for a given process when the number of final particles increases.

For example, we have only 3 diagrams for e+e- --> W+ W- , but when one photon is added, e+e- --> W+ W- gamma , then 18 diagrams appear even after omitting the tiny interaction between e and scalar bosons( Higgs and Goldstone bosons ). Addition of one another photon, e+e- --> W+ W- gamma gamma , yields 138 diagrams. Further if one wants to make more realistic calculation around the threshold of W pair production, taking into account the decay of W , say, W- --> e- nu_e-bar and W+ --> u-dbar , then one has to consider 24 diagrams for e+e- --> e- nu_e-bar u d-bar and 202 for e+e- --> e- nu_e-bar u d-bar gamma . In unitary gauge, as only physical particles appear in the Lagrangian, the numbers of diagrams are less than those mentioned above.

One may think that it is enough to select several diagrams which dominate the cross section. Even if one can find such dominant diagrams, one has to respect the gauge invariance among this subset of diagrams. Usually number of diagrams in the gauge invariant subset is not so small. For example, for the process e+e- --> nu_e nu_e-bar W+ W- , we have 60 diagrams in all. Among them 30 diagrams form one gauge invariant set and the rest does another one. Hence still we meet the same difficulty to handle with many diagrams. In addition, there remains a possibility that the experimental cuts imposed on the final particles renders the dominant diagrams to be less prominent and all diagrams give somehow the same order of magnitude to the cross section. If this is the case, one has to keep whole the diagrams in the calculation after all.

Through the numerous experiments done at e+e- colliders, we have learned that higher order corrections should be included when we want to compare theories with experimental data in detail. This implies that we have to calculate at least one-loop corrections to a given process. As an example, consider the process e+e- --> W+ W- gamma . To regularize the infrared divergence due to soft photon emission, we have to include loop diagrams for e+e- --> W+ W- beyond the tree level, which contain virtual photon exchange and remove the divergence when combined with real photon emission process. The requirement of gauge invariance among one-loop diagrams demands, in turn, inclusion of other one-loop diagrams with exchanges of Z0 , W or other possible particles. Then it is clear that the total number of diagrams becomes very huge and it is almost impossible even to enumerate all diagrams. In many cases it seems out of ability of mankind. For simple W-pair production, in general covariant gauge, the number is around 200 diagrams in the same approximation stated above, but for e- nu_e-bar u d-bar it amounts more than 3,700.

Facing to the difficulty described above, we cannot help to find some ways to get rid of. As a solution we can choose the following one: As diagrams are constructed based on a set of definite rules , Feynman rules, it is natural to develop a computer code which can generate all the diagrams to any process, once initial and final particles are given. It should be able not only to enumerate diagrams but also generate automatically relevant amplitudes to be evaluated on computers, in other words, create a FORTRAN source code ready for amplitude calculation. GRACE is such a system that realizes this idea and help us to reduce the most tedious part of works.

- Specify the process.
- Choose appropriate models.
- Fix the order of perturbation( at the tree level, this is unique ).
- Enumerate all possible diagrams.
- Write down amplitudes.
- Prepare the kinematics for final particles.
- Integrate the amplitude squared in the phase space of final particles, including experimental cuts, if necessary.
- Generate events so that the simulation of the process in a detector is available.
- Check the results.

- All the tree diagrams for a given process when the order of perturbation is fixed.
- Diagrams on X-window and its print-out.
- FORTRAN source code which contains helicity amplitude of the process.
- Default values of all physical constants, except for the strong coupling constant .
- Interface routines to the program package CHANEL which contains subroutines designed to evaluate the amplitude.
- Default code for kinematics.
- Interface routines to the multi-dimensional integration package BASES.
- Interface routines to the event generation package SPRING.
- Test program for gauge invariance check of the generated amplitude.
- Any diagram and its amplitude can be omitted in the calculation by setting the appropriate flags off. In the integration step the unitary gauge is the default.

- names of initial particles,
- names of final particles,
- order of perturbation in QED, electroweak and/or QCD.
- code number of kinematics.

When GARCE is initiated with the data file containing these inputs, it constructs all possible diagrams and creates an output file to draw all Feynman diagrams for the convenience of the user to look them by eyes. At the same time a set of FORTRAN subprograms is generated. These include those which are needed to calculate the amplitude with the help of CHANEL , to integrate over phase space by BASES and to generate events by SPRING .

After all the programs are successfully generated, the remaining tasks for user are

- to examine the kinematics,
- to check and to edit some parameters in a few subroutines,
- to check the gauge invariance of the amplitude,
- to supply the value of strong coupling constant.

- Graph generation subsystem

When initial and final states of the elementary process are given as the input as well as the orders of couplings, a complete set of Feynman graphs is generated according to the theoretical model defined in a model definition file. For the time being QED, Electroweak and QCD models in the tree level are supported. The information of generated graphs is stored in a file as an output.

Reading the graph information from the file, the graph drawer displays the Feynman graphs on the screen under the X-Window system or prints them on a paper. - Source generation subsystem

From the graph information produced by the first subsystem, a FORTRAN source code is generated in a form of program components suited for the numerical integration package BASES and the event generation package SPRING .

The source code is constructed based on our helicity amplitude formalism, which consists of many calling sequences of subprograms given in CHANEL and its interface routines.

In addition to these program components, the subsystem generates a main program, by which the gauge invariance of the generated amplitudes can be tested. - Numerical integration subsystem

Combining the generated source code together with the kinematics routines and the GRACE library, the numerical integration is performed by BASES to obtain the total cross section. For this, however, in general, one may have to prepare the kinematics routines when the default one is not appropriate. As the output of integration, the numerical value of total cross section, the convergency behavior of integration, one and two dimensional distributions of the cross section are given besides the probability information in a file, which is used in the event generation. Looking the convergency behavior carefully one can judge if the resultant value is reliable or not. - Event generation subsystem

Using almost all the same subprograms in the integration, events with weight one are generated by the event generation program SPRING . To achieve a high generation efficiency, it uses the probability information produced by BASES . Conceptually, SPRING samples a point in the integration volume according to the probability. If the probability information is a complete one, the sampled point is exactly corresponding to a generating event. Since, however, it is impossible to get a complete information numerically, the sampled point is tested whether it is accepted or not.

The user can record and analyze the generated events according to the information in the following chapters.

X(I), I = 1, .... , NDIM ,

and assume their values are normalized in, say, [0,1] . Then we have to translate these variables into four-momentum of final particle, say J -th particle, P(1,J), P(2,J), P(3,J), P(4,J) of total N particles ( in GRACE , P(4,J) is the energy ),

X(I) --> P(K,J) . K =1,...,4, J=1,..., N

This is known as kinematics for the given process. This mapping is not always unique and in some cases a single value of X(I) may correspond to multi-value of particle momenta. GRACE , unfortunately, only can give a candidate of the kinematics. The reason is that the present popular integration packages, such as BASES or VEGAS , utilize a special algorithm to search for the singularities of the integrand. The matrix element squared, the integrand, becomes singular when the denominators of propagators of internal particles become very small compared with the typical energy of the process considered. This happens when a mass of an internal line is very small. As is well known, if a singularity is running along the diagonal in a plane of two integration variables, these programs cannot give reliable estimate of the integral, because they fail to catch the singularity at all. In order to get good convergence of the integration over many iterations, all the singularities must be parallel to the integration axes. This means that these peaks located in the space of kinematical variables, are mapped onto the line of constant value of some X(I) . In order to do this, we have to choose very carefully the transformation between random numbers and kinematical variables. The typical kinds of singularities we meet in real calculation are as follows;

- mass singularity
- infrared singularity
- t -channel photon exchange
- resonance formation( decay of heavy particles )

In some processes the number of independent variables is greater than that of singularities, and one can easily find a kinematics which is suitable to make them smooth. If this is not the case, however, one may not be able to find such good kinematics to avoid diagonal singularity even after much efforts. Hence it is quite difficult to give the general kinematics which is capable of dealing with all kinds of singularities at once, or a single set of transformations.

Considering the situation described above, GRACE only provides some sample of kinematics to users. So the following points are left for the users:

- Select the ordering of particles to use kinematics.
- Select the kinematics among the candidates.
- Set physical parameters, W=center-of-mass energy.
- Change the parameters for the integration, e.g., the number of iteration, number of required accuracy and so on.
- Revise or write the source code by the user when all candidates provided by the system are no good.