ALHEP symbolic algebra program for high-energy physics
ALHEP symbolic algebra program for
high-energy physics
V. Makarenko 1
NC PHEP BSU, 153 Bogdanovicha str., 220040 Minsk, Belarus
Abstract
ALHEP is the symbolic algebra program for high-energy physics. It deals with
amplitudes calculation, matrix element squaring, Wick theorem, dimensional regu-
larization, tensor reduction of loop integrals and simplification of final expressions.
The program output includes: Fortran code for differential cross section, Mathe-
matica files to view results and intermediate steps and TeX source for Feynman
diagrams. The PYTHIA interface is available.
The project website http://www.hep.by/alhep contains up-to-date executa-
bles, manual and script examples.
1 Introduction
The analytical calculations in high-energy physics are mostly impossible with-
out a powerful computing tool. The big variety of packages is commonly used
[1,2,3,4,5,6,7,8,9,10,11,12,13,14]. Some are general-purpose symbolic algebra
programs with specific HEP-related plug-ins (REDUCE [15], Mathematica
[16]), some are designed especially for particle physics (CompHEP [1], SANC
[2], GRACE [3] etc.) and some are created for specific interaction class or spe-
cific task. Many of them uses external symbolic algebra core (Form [17], Math-
Link [16]). They can deal with matrix elements squaring (FeynCalc [7]) or cal-
culate helicity amplitudes directly (MadGraph [4], CompHEP [1], O’Mega [5]).
Some packages provide numerical calculations, some require external Monte-
Carlo generator to be linked. Some programs contain also one-loop calculation
routines (FormCalc [6], GRACE [3], SANC [2]). Nevertheless, there is no uni-
form program that meets all the user requirements.
Email address: makarenko@hep.by (V. Makarenko ).
1 Supported by INTAS YS Grant no. 05-112-5429
Preprint submitted to Elsevier 17 June 2021
http://arxiv.org/abs/0704.1839v1
Every calculation requires the program-independent check. The optimal tac-
tics is the simultaneous usage of two (or more) different symbolic algebra
packages.
ALHEP is a symbolic algebra program for performing the way from Standard
Model Lagrangian to amplitude or squared matrix element for the specified
scattering process. It can also be useful for loop diagrams analysis. The basic
features are:
• Diagrams generation using Wick theorem and SM Lagrangian.
• Amplitude calculation or matrix element squaring.
• Bondarev functions method for traces calculation.
• Tensor reduction of loop integrals.
• Dimensional regularization scheme.
• Generation of Fortran procedures for numerical analysis (PYTHIA [18] and
LoopTools [19] interfaces are implemented).
The current ALHEP version have several implementation restrictions, that will
be lifted in future. The following features are in progress of implementation:
• Bremsstrahlung part of radiative correction (the integration over real photon
phase space).
• Complete one-loop renormalization scheme including renormalization con-
stants derivation.
• Arbitrary Lagrangian assignment.
After these methods implementation the complete one-loop analysis will be
available. Please refer to project website for program updates.
ALHEP website http://www.hep.by/alhep contains the up-to-date executa-
bles (for both Linux & Win32 platforms), manual and script examples. The
mirror at http://cern.ch/~makarenko/alhep is also updated.
2 ALHEP Review. Program Structure
The ALHEP program internal structure can be outlined as follows:
• The native symbolic algebra core.
• Common algebra libraries:
Dirac matrices, tensor and spinor algebra, field operators & particle wave
functions zoo.
• Specific HEP functions and libraries.
It include Feynman diagrams generation, trace calculations, helicity am-
plitudes method, HEP-specific simplification procedures, tensor integrals
reduction and others.
• Interfaces to Mathematica, Fortran, TeX and internal I/O format.
Fortran code is used for further numerical analysis. Mathematica code can
be used for viewing any symbolic expression in program. But no backward
interface from Mathematica is currently implemented. TeX output can be
generated for Feynman diagrams view. Internal I/O format is implemented
for the most of symbolic expressions, allowing save and restore calculation
at intermediate steps.
• Command script processor.
User interface is implemented in terms of command scripts. The ALHEP
script language have C-like syntax, variables, arithmetic operations, func-
tion calls. All HEP-related tasks are implemented as build-in functions.
2.1 Getting Started
To use ALHEP one should download the pre-compiled executable for appro-
priate platform (Linux and Win32 are available) and write a control script
to describe your task. ALHEP program should be launched with the single
argument: script file name to be invoked The following steps are required to
create a workspace:
• Download ALHEP executables from project website: alhep.gz (for Linux)
or alhep.zip (Win32). Unpack executable, e.g.
gzip -d alhep.gz
• Download up-to date command list: ALHEPCommands.txt. The set of com-
mands (or options) may be changed in future versions and this manual may
be somewhat obsolete. Please refer to ALHEPCommands.txt file that always
corresponds to the latest ALHEP version. The last changes are outlined in
RecentChanges.txt file at website.
• Create some working directory and compose command script file therein.
For example consider the uū → W+W−γ process with {−+}-helicities of
initial quarks. To calculate amplitude we create the following script file (call
it "test.al"):
SetKinematics(2, 3 // 2->3 process
,QUARK_U,"p\_1","e\_1" // u
,-QUARK_U,"p\_2","e\_2" // u-bar
,WBOZON,"f\_1","g\_1" // W{+}
,PHOTON,"f\_0","g\_0" // photon
,-WBOZON, "f\_2", "g\_2" ); // W{-}
SetDiagramPhysics(PHYS_SM_Q1GEN);// SM with 2 quarks only
SetMassCalcOrder(QUARK_U, 0); // consider massless
SetMassCalcOrder(QUARK_D, 0); // consider massless
diags = ComposeDiagrams(3); // create diagrams, e^3
DrawDiagrams(diags, "res.tex",
DD_SMALL|DD_SWAP_TALES, FILE_START);
SetFermionHelicity(1, -1);
SetFermionHelicity(2, 1);
SetParameter(PAR_TRACES_BONDEREV, 1);
ampl = CalcAmplitude(RetrieveME(diags));
ampl = KinArrange(ampl); // arrange result
ampl = Minimize(ampl); // minimize result
SaveNB("res.nb",ampl,"",FILE_START|FILE_CLOSE);
f = NewFortranFile("res.F", CODE_F77);
CreateFortranProc(f, "UUWWA", ampl,
CODE_IS_AMPLITUDE|CODE_COMPLEX16
|CODE_CHECK_DENOMINATORS|CODE_PYTHIA_VECTORS);
The SetKinematics() function declares particles, momenta and polariza-
tion symbols. The physics is declared with PHYS_SM_Q1GEN option to restrict
diagrams number. The amplitude for all (d, c, b) internal quarks can be ob-
tained from generated here by quark mixing matrix replacement: U2ud →
U2ud+U
ub (chiral limit). Diagrams are created by ComposeDiagrams()
call. CalcAmplitude() function creates the symbolic value for process am-
plitude. For detailed discussion of this example see sec. 15.1.
• Create simple batch file "run.me" like:
~/alhep_bin_path/alhep test.al
The ALHEP program creates some console output (current commands,
scroll bars and some debugging data). If it is not allowed one should redirect
console output to file here.
The batch execution with test.al command file takes about 1 minute at
1.8GHz P4 processor. The following files are created:
res.nb: Mathematica file containing symbolic expression of amplitude. Cre-
ated by SaveNB() function.
res.F: F77 code for numerical analysis created by NewFortranFile() and
CreateFortranProc() calls. The library file alhep_lib.F is required for
code compilation and should be downloaded from project website. See sec.
11.1 for Fortran generation and compilation review.
res.tex: TeX source for 19 Feynman graphs generated. The AxoDraw [23]
Latex package is used. The res.tex file should be included into your La-
TeX document using \input res.tex command. The template document
to include your diagrams can be found at program website. See fig. 1 in sec.
15.1 for diagrams generated.
debug.nb: Mathematica file with debugging information and some interme-
diate steps. The amount of debugging information is declared in debug.ini
file in working directory. See sec. 13 for further details.
See sec. 15 and project website for another examples. It is convenient to use
some example script as template and modify it for your purposes.
2.2 Calculation scheme
The usual ALHEP script contains several steps:
(1) Initialization section
Declaration of process kinematics: initial and final-state particles, titles
for particle momenta & polarization vectors (SetKinematics()).
Physics model definition. The SM physics or part of SM Hamiltonian
should be specified (SetDiagramPhysics()). The shorter Hamiltonian is
selected, the faster is Wick theorem invocation.
Polarization declaration. Every particle is considered as polarized with
abstract polarization vector by default. The specific helicity value can be
set manually or particle can be marked as unpolarized. The several ways
of polarization involving are inplemented, see sec. 4.4 (SetPolarized(),
SetFermionHelicity(), SetPhotonHelicity(), ...).
Setting mass-order rules for specific particles (SetMassCalcOrder()).
One can demand the massless calculation for light particle, that greatly
saves evaluation time. One can also demand keeping particle mass with
specific order Mn and drop out the higher-order expressions like Mn+1.
It allows to consider the leading mass contribution without calculating
precisely.
(2) Diagrams generation
Feynman diagrams are generated for specific en order using the Wick
theorem algorithm (ComposeDiagrams()) User may draw diagrams here
(DrawDiagrams()), halt the program (Halt()) and check out if diagrams
are generated correctly.
After the diagram set is generated one may cut-off not interesting di-
agrams to work with the shorter set or select the single diagram to work
with (SelectDiagrams()). The loop corrections are calculating faster
when processed by single diagrams.
Then matrix element is retrieved from diagrams set (RetrieveME()).
Before any operation with loop matrix element one should declare
the N−dimensional space (SetNDimensionSpace()). Some procedures
involve N−dimensional mode automatically for loop objects, but most
functions (arranging, simplification etc.) don’t know the nature of expres-
sion they work with. Therefore the dimensional mode should be forced.
The diagrams set and any symbolic expression may be saved and re-
stored at next session to refrain from job repetition (Save(), Load()).
(3) Amplitude calculation
The amplitude evaluation (CalcAmplitude()) is the faster way for
multi-diagram process analysis.
All the particles are considered polarized. The spinor objects are pro-
jected to abstract basis spinor. The basis vectors are generated in nu-
merical code to meet the non-zero denominators condition. See sec. 7 for
method details.
(4) Matrix element squaring (coupling to other)
The squaring procedure (SquareME()) is controlled by plenty of op-
tions, intended mostly for the performance tuning and debugging.
It basically includes reduction of gamma-matrix sequences coupled to
kinematically dependent vectors. It reduces the number of matrices in
every product to minimum. The item-by-item squaring is followed. For
loop × born∗ couplings the virtual integrals are involved. See sec. 8 for
details.
Amplitude calculation is fast, but its result may be more complicated
than squaring expressions. Amplitude depends on particle momenta, po-
larization vectors and additional basis vectors. For unpolarized process
the averaging cycle is generated in Fortran code, and complex-numbers
calculation should be performed. The squaring result for unpolarized pro-
cess is the polynomial of momenta couplings only. Hence there is no
unique answer what result is simpler for 10-15 diagrams reaction. One
should definitely use amplitude method for more than 10-15 diagrams
squaring.
(5) Loop diagrams analysis
The tensor virtual integrals are reduced to scalar ones (Evaluate(),
sec. 9). The scalar coefficients of tensor integral decompositions are also
reduced to scalar integrals in the most cases (for 1− 4 point integrals).
For scalar loop integrals the tabulated values are used. There is no
reason to tabulate integrals with complicated structure. Hence scalar in-
tegrals table contain the A0 and B0 integrals with different mass configu-
ration. It also contains a useful D0 chiral decomposition. Other integrals
should to be resolved using LoopTools-like [19,20] numerical programs.
The renormalization procedure is under construction now. The counter-
terms (CT) part of Lagrangian leads to CT diagrams generating. In the
nearest future the abstract renormalization constants (δm, δf etc.) will be
involved and tabulated for the minimal on-shell scheme. The automatic
derivation of constants is supposed to be implemented further. Please
refer ALHEP website for implementation progress.
(6) Simplification
The kinematic simplification procedure is available (KinSimplify(),
sec. 12). It reduces expression using all the possible kinematic relations
between momenta and invariants. The minimization of +/* operations
in huge expressions can also be performed (Minimize()).
(7) Fortran procedure creation for numerical analysis
F90 or F77-syntaxes for generated procedures are used. Generated code
can be linked to PYTHIA, LoopTools and any Monte-Carlo generator for
numerical analysis.
3 ALHEP script language
The script syntax is similar to C/Java languages. Command line breaks at
”;” symbol only. Comments are marked as // or /*...*/. Operands may be
variables or function calls. If no function with some title is defined, the operand
is considered as variable. The notation is case-sensitive.
The script language have no user-defined functions, classes or loop operators.
It seems to be useless in current version.
All ALHEP features are implemented as build-in functions. The execution
starts from the first line and finishes at the end of file (or at Halt() command).
Variable types are casted and checked automatically (in run-time mode). No
manual type specifying or casting are available. List of ALHEP script internal
types:
• Abstract Symbolic Expression (expr).
Result of function operations. Can be stored to file(Save()) and loaded
back(Load()). Basic operations: +,-,*. The division is implemented using
Frac(a,b) function. Are not supposed to be inputed manually, although a
few commands for manual input exist.
• Integer Value (int).
Any number parameter will be casted to integer value. Basic operations:
+,-,*,|. The division is performed using Frac(a,b) function. No fractional
values are currently available.
• String (str).
String parameters are started and closed by double quotes (”string vari-
able”). Used to specify symbol notation (e.g. momenta titles), file names,
etc. Basic operations: +(concatenate strings).
• Set of Feynman Diagrams (diagrams).
Result of diagrams composing function. Basic operations: Save(), Load()
and SelectDiagrams(). Diagrams can also be converted into TeX graphics.
• Matrix Element (me).
Contains symbolic expression for matrix element and information on it’s
use: list of virtual momenta (integration list) etc.
The total up-to-date functions list can be found in ALHEPCommands.txt file
at project website.
4 Initialization section
4.1 Particles
Particles are determined by integer number, called particle kind (PK). The
following integer constants are defined: ELECTRON, MUON, TAULEPTON, PHOTON,
ZBOZON, WBOZON, QUARK_D, QUARK_U, QUARK_S, QUARK_C, QUARK_B, QUARK_T,
NEUTRINO_ELECTRON, NEUTRINO_MUON, NEUTRINO_TAU. The ghost and scalar
particles are not supposed to be external and their codes are unavailable.
Antiparticles have negative PK that is obtained by ”-PK” operation.
If kinematics is declared the particles can be secelted using the particle ID
(PID) number. Initial particles have negeative ID (-1, -2) and final are pointed
by positive numbers (1, 2, 3...).
4.2 Kinematic selection
Before any computations may be performed one needs to declare the kinematic
conditions. They are: number of initial and final state particles, PK codes,
symbols for momentum and polarization vector for every particle.
SetKinematics((int)N Initials, (int)N Finals,
[(int)PK I, (str)momentum I, (str)polarization I, ...]);
Here N_Initials and N_Finals – numbers of initial and final particles in kine-
matics. The next parameters are particle kind, momentum and polarization
symbols, repeated for every particle. For example, the e−(k1, e1)e
+(k2, e2) →
µ−(p1, e3)µ
+(p2, e4) process should be declared as follows:
SetKinematics(2, 2,
ELECTRON,"k\_1","e\_1",
-ELECTRON,"k\_2","e\_2",
MUON, "p\_1", "e\_3",
-MUON, "p\_2", "e\_4");
4.3 Particle masses
All the particles are considered massive by default. Particle can be declared
massless using the SetMassCalcOrder function.
SetMassCalcOrder((int)PK, (int)order);
PK: particle kind,
order: maximum order of particle mass to be kept in calculations. Zero value
means massless calculations for specified particle. Negative value declares
the mass-exact operations.
The mass symbols are generated automatically and look like me, mW etc.
Hence, all the electrons in process will have the same mass symbol. To declare
unique mass symbols for specific particles the SetMassSym() function is used.
The different masses for unique particles are often required to involve the
Breit-Wigner distribution for particles masses.
SetMassSym((int)PID, (str or expr)mass);
PID: particle ID in kinematics (< 0 for initial and > 0 for final particles),
mass: new mass symbol, like "m\_X" for mX .
4.4 Polarization data
All the particles are considered as polarized initially. The default polarization
vector symbols are set together with kinematic data (in SetKinematics()).
To declare particle unpolarized the SetPolarized() function is used:
SetPolarized((int)PID, (int)polarized);
PID: particle ID in kinematics,
polarized: 1 (polarized) or 0 (unpolarized).
One can declare the specific polarization state (helicity value) for particles.
The photon (k,e) helicities are involved in terms of two outer physical mo-
menta:
e±µ =
(p+ · k)p−µ − (p− · k)p+µ ± iǫµαβνpα+p
/Nk,p±. (1)
SetPhotonHelicity((int)PID, (int)h, (str)”pP”, (str)”pM”]);
PID: photon particle ID in kinematics,
h: helicity value: ±1 or 0. Zero value clears the helicity information
pP, pM: p± base vectors in (1) formula. The γ-coupled forms of (1) (precise
or chiral [26]) are used automatically if available.
One may select the transverse unpolarized photon density matrix (instead of
usual −(1/2)gµν):
ν → −gµν +
(P · k)
(P · k)
kµkνk
(P · k)2
. (2)
SetPhotonDMBase((int)PID, (str)”P”);
PID: photon particle ID in kinematics,
"P": Basis vector P in (2).
The fermion (k,e) helicity is declared using:
SetFermionHelicity((int)PID, (int)h);
SetFermionHelicity((int)PID, (srt)”h”);
PID: fermion particle ID in kinematics,
h: helicity value: ±1 or 0. Zero value clears the helicity information. Density
matrix is usual:
uū → p̂γ± (massless fermion).
"h": symbol for scalar parameter in the following density matrix (massless
fermion):
uū → (1/2)p̂(1 + hγ5).
Notes:
• SetPolarized(PID, 0) call will also clear helicity data.
• SetXXXHelicity(PID, 1, ...) also sets particle as polarized (previous
SetPolarized(PID,0) call is canceled).
• SetXXXHelicity(PID, 0, ...) clears helicity data but does not set parti-
cle unpolarized.
4.5 Physics selection
One can either declare the full Standard Model (in unitary or Feynman gauge)
or select the parts of SM Lagrangian to be used. The QED physics is used by
default. The Feynman rules corresponds to [27] paper.
SetDiagramPhysics( (int)physID );
physID: physics descriptor, to be constructed from the following flags:
PHYS_QED: Pure QED interactions,
PHYS_Z and PHYS_W: Z- and W-boson vertices,
PHYS_MU, PHYS_TAU: Muons and tau leptons (and corresponding neutrinos),
PHYS_SCALARS: Scalar particles including Higgs bosons,
PHYS_GHOSTS: Faddeev-Popov ghosts,
PHYS_GAUGE_UNITARY: Use unitary gauge,
PHYS_QUARKS, PHYS_QUARKS_2GEN and PHYS_QUARKS_3GEN: {d, u}, {d, u, s, c}
and {d, u, s, c, b, t} sets of quarks,
PHYS_RARE_VERICES: vertices with 3 or more SCALAR/GHOST tales,
PHYS_CT: Renormalization counter-terms (implementation in progress),
PHYS_SM: full SM physics in unitary gauge (all the flags above except for
PHYS_CT),
PHYS_SM_Q1GEN: the SM physics in unitary gauge with only first generation
of quarks (d, u),
PHYS_ELW: SM in Feynman gauge with no rare 3-scalar vertice,
PHYS_EONLY: no muons, tau-leptons and adjacent neutrinos,
PHYS_NOQUARKS: no quarks,
PHYS_4BOSONS_ANOMALOUS: anomalous quartic gauge boson interactions (see
[28]), affects on AAWW and AZWW vertices.
The less items are selected in Lagrangian, the faster diagram generation is
performed.
5 Bondarev functions
The Bondarev method of trace calculation is implemented according to [25]
paper. The trace of γµ-matrices product becomes much shorter in terms of
F -functions. The number of items for Tr[(1 − γ5)â1â2 · · · â2n] occurs 2n. For
example, the 12 matrices trace contains 10395 items usually while the new
method leads to 64-items sum.
The 8 complex functions are introduced:
F1(a, b) = 2[(aq−)(bq+)− (ae+)(be−)],
F2(a, b) = 2[(ae+)(bq−)− (aq−)(be+)],
F3(a, b) = 2[(aq+)(bq−)− (ae−)(be+)],
F4(a, b) = 2[(ae−)(bq+)− (aq+)(be−)], (3)
F5(a, b) = 2[(aq−)(bq+)− (ae−)(be+)],
F6(a, b) = 2[(ae−)(bq−)− (aq−)(be−)],
F7(a, b) = 2[(aq+)(bq−)− (ae+)(be−)],
F8(a, b) = 2[(ae+)(bq+)− (aq+)(be+)].
The basis vectors q± and e± are selected as follows:
(1, ±1, 0, 0), eµ± =
(0, 0, 1,±i). (4)
The results for traces evaluation looks as follows:
Tr[(1− γ5)â1â2] = F1(a1, a2) + F3(a1, a2),
T r[(1− γ5)â1â2â3â4] = F1(a1, a2)F1(a3, a4) + F2(a1, a2)F4(a3, a4) +
+F3(a1, a2)F3(a3, a4) + F4(a1, a2)F2(a3, a4).
Please refer to [25] paper for method details.
SetParameter(PAR TRACES BONDEREV, (int)par);
par: 1 or 0 – allow or forbid Bondarev functions usage.
The numerical code for F -functions (3) is contained in alhep_lib.F library
file. The code is available for scalar couplings only: Fµνp
µqν . If vector Bondarev
functions remain in result, the Fortran-generation procedure fails. One should
repeat the whole calculation without Bondarev functions in that case.
6 Diagrams generation
The diagrams are generated after the kinematics and physics model are de-
clared. The Wick-theorem-based method is implemented.
The only distinct from the usual Feynman diagrams is the following: the vertex
rules have additional 1/2 factors for every identical lines pair. It makes two
effects:
- The crossing diagrams are usually involved if they have different topology
from original. I.e. if two external photon lines starts from single vertex, they
are not crossed. Nevertheless all the similar external lines have crossings in
ALHEP.
- Some ALHEP diagrams have 2n factors due to identical internal lines. For
example, the W+W− → {γγWWvertex} → γγ → {γγWWvertex} →
W+W− diagram will have additional factor 4. The result remains correct
due to 1/2 factors at every γγW+W− vertex.
The diagrams are generated without any crossings. The crossed graphs are
added automatically during the squaring or amplitude calculation procedures.
diagrams ComposeDiagrams((int)n);
n: order of diagrams to be created, MX→Y ∼ en.
Uses current physics and kinematics information. Returns generated diagrams
One may select specified diagrams into another diagrams set:
diagrams SelectDiagrams( (diagrams)d, (int)i0 [, i1, i2...]);
d: initial diagrams set. Remains unaffected during the procedure.
i0..iN: numbers of diagrams to be selected. First diagram is ”0”.
The new diagrams set is returned.
To retrieve matrix element from the diagram set the RetrieveME() is used.
me RetrieveME( (diagrams)d );
d: diagrams list.
7 Helicity amplitudes
The amplitudes are calculated according to [24] paper. Every spinor in matrix
element is projected to common abstract spinor:
ūp =
ūQupūp
ūQup
= eiC
ūQPp
(Tr[PpPQ])1/2
, up =
upūpuQ
ūpuQ
= eiC
(Tr[PQPp])1/2
The eiC factor is equal for all diagrams and may be neglected.
The projection operator PQ = upūp is choosen as follows:
• PQ = Q̂(1 + γ5)/2 – for massive external fermions,
• PQ = Q̂(1 + ÊQ)/2 – if one of fermions is massless.
The value of vector Q is selected arbitrary in Fortran numerical procedure.
The additional basis vector EQ (if exists) is selected meet the polarization
requirements ((EQ.Q) = 0, (EQ.EQ) = −1). The fractions like 1/Tr[PQPp]
may turn 1/0 at some Q and EQ values. The denominators check procedures
are generated in Fortran code (the CODE_CHECK_DENOMINATORS key should be
used in CreateFortranProc() call). If the | Tr[PQPp] |> δ check is failed, the
another Q and EQ values are generated.
expr CalcAmplitude( (me)ME );
ME: matrix element retrieved from the whole diagrams set.
The result expression is a function of all the particle helicity vectors. The
averaging over polarization vectors is performed numerically. The numerical
averaging procedure is automatically generated in Fortran output if unpolar-
izaed particles are declared.
The CODE_IS_AMPLITUDE key in CreateFortranProc() procedure declares
expression as amplitude and leads to proper numerical code (Ampl×Amlp∗).
8 Matrix Element squaring
The squaring procedure have the following steps:
• matrix elements simplification to minimize the γ-matrices number,
• denominators caching to make procedure faster,
• item-by-item squaring (coupling to other conjugated),
• saving memory mechanism to avoid huge sums arranging (SQR_SAVE_MEMORY
option),
• virtual integrals reconstruction.
expr SquareME((me)ME1, [(me)ME2,] [(int)flags ]);
ME1: matrix element #1,
ME2: matrix element #2 (should be omitted for squaring),
flags: method options (defauls is ”0”):
SQR_CMS: c.m.s. consideration (initial momenta are collinear). The additional
pseudo-covariant relations (p1.ε2 → 0, p2.ε1 → 0) appear that simplify work
with abstract polarization vectors.
SQR_NO_CROSSING_1: do not involve crossings for ME #1,
SQR_NO_CROSSING_2: do not involve crossings for ME #2,
SQR_MANDELSTAMS: allow Mandelstam variable usage,
SQR_PH_GAMMA_CHIRAL: tries to involve photons helicities in short chiral ê±
form (see sec. 4.4, [26]).
SQR_PH_GAMMA_PRECISE: involve precise ê± form for photon helicities. If the
polarization vectors are not coupled to γµ, the vector form e
µ is used (see
sec. 4.4).
SQR_SAVE_MEMORY: save memory and processor time for huge matrix ele-
ment squaring. Minimizes sub-results by every 1000 items and skips final
arranging of the whole sum. No huge sums occurs in calculation in this
mode, but the result is also not minimal. If result expression contains a
sum of 105 − 106 items (when expanded), the arranging time is significant
and SQR_SAVE_MEMORY flag should be involved. If no results are calculated
in reasonable time the CalcAmplitude() (see sec. 7) procedure should be
used.
If two matrix elements are given, the first will be conjugated, i.e. the result is
ME1∗ ×ME2.
9 Virtual integrals operations
The tensor virtual integrals are reduced to scalar ones using two methods.
If tensor integral is coupled to external momentum Iµp
µ and p-vector can
be decomposed by integral vector parameters, the fast reduction is involved.
The Dx integrals for 2 → 2 process contain the whole basis of 4-dimension
space and Dx-couplings to any external momentum can be decomposed. It
works well if all the polarization vectors are constructed in terms of external
momenta.
The common tensor reduction scheme is involved elsewhere. Tensor integrals
are decomposed by the vector basis like
Iµν (p, q) = I
pµpν , qνqµ, p(µqν), gµν
The linear system for scalar coefficients is composed and solved. Implemented
forBj and Ci,ij integrals only. The other scalar coefficients should be calculated
numerically [19].
expr Evaluate( (expr)src );
Evaluate: Reduction of tensor virtual integrals to scalar ones. The scalar
coefficients in tensor VI decomposition are also evaluated (not for all integrals).
expr ConvertInvariantVI( (expr)src );
ConvertInvariantVI: Vector parameters are substituted by scalars. I.e.:
C0(k1, k2, m1, m2, m3) → C0(k21, (k1 − k2)2, k22, m1, m2, m3).
Note: A- and B- integrals are converted automatically during arrangement.
expr CalcScalarVI( (expr)src );
CalcScalarVI: Substitutes known scalar integrals with its values. The most
of UV-divergent integrals (A0,B0) are substituted. The chiral decomposition
for D0-integral is also applied. The complicated integrals should be calculated
numerically [19].
The source expressions are unaffected in all the functions above.
10 Regularization
The dimensional regularization scheme is implemented. One may change the
space-time dimension before every operation:
SetNDimensionSpace( (int)val );
val: 1 (n-dimensions) or 0 (4-dimensions space).
expr SingularArrange( (expr)src );
SingularArrange: turns expression to 4-dimensional form. Calculates (n−4)i
factor in every item and drops out all the neglecting contributions.
SetDRMassPK( (int)PK )
SetDRMassPK: set the particle to be used as DR mass regulator.
PK: particle kind (see sec. 4.1). ”0” value declares the default ”µ” DR mass.
11 ALHEP interfaces
11.1 Fortran numerical code
The numerical analysis in particle physics is commonly performed using the
Fortran programming language. Hence, we should provide the Fortran code to
meet the variety of existing Monte-Carlo generators.
To start a new Fortran file the NewFortranFile() function is used:
int NewFortranFile( (str)fn [, (int)type]);
fn: output file name
type: Fortran compiler conventions: CODE_F77 or CODE_F90. The FORTRAN 77
conventions are presumed by default.
Function returns the ID of created file.
Then we can add a function to FORTAN file:
CreateFortranProc( (int)fID, (str)name, (expr)src [, (int)keys]);
fID: file ID returned by NewFortranFile() call.
name: Fortran function name.
src: symbolic expression to be calculated. The source may be | M |2, dσ/dΓ
(use CODE_IS_DIFF_CS flag) or amplitude (CODE_IS_AMPLITUDE flag). The
result is always dσ/dΓ calculation procedure.
keys: option flags for code generation (default is CODE_REAL8):
CODE_REAL8: mean all the symbols in expression as REAL*8 values.
CODE_COMPLEX16: declare variables type as COMPLEX*16.
CODE_IS_DIFF_CS: the source expression is differential cross section.
CODE_IS_AMPLITUDE the source expression is amplitude and a squaring code
should be generated: AMPL*DCONJG(AMPL).
CODE_CHECK_DENOMINATORS check denominators for zero. Used to re-generate
free basis vectors in amplitude code (see sec. 7).
CODE_LOOPTOOLS: Create LoopTools [19] calls for virtual integrals.
CODE_SEPARATE_VI create unique title for every virtual integral function (ac-
cording to parameter values). Do not use with CODE_LOOPTOOLS.
CODE_PYTHIA_VECTORS retrieve vector values from PYTHIA [18] user process
PUP(I,J) array.
CODE_POWER_PARAMS factorize and pre-calculate powers if possible.
CODE_NO_4VEC_ENTRY do not create a 4-vector entry for function.
CODE_NO_CONSTANTS do not use predefined physics constants. All variables
becomes external parameters.
CODE_NO_SPLIT do not split functions by 100 lines.
CODE_NO_COMMON: don’t use CONNON-block to keep internal variables.
The scalar vector couplings, Bondarev functions (see sec. 5) and εabcdp
objects are replaced with scalar parameters to be calculated once. These func-
tions are calculated using alhep_lib.F library procedures.
Some compilers works extremely slow with long procedures. Therefore Fortran
functions are automatically splitted after every 100 lines for faster compilation.
The internal functions functions have TEMPXXX() notation. To avoid problems
with several ALHEP-genarated files linking we should rename the TEMP-prefix
to make internal functions unique.
SetParameterS(PAR STR FORTRAN TEMP, (str)prefix);
prefix: "TMP1" or another unique prefix for temporary functions name.
To obtain the better performance of numerical calculations ALHEP provides
the mechanism for minimization of +/× operations in the expression. We
recommend to invoke the Minimize() function before Fortran generation (see
sec. 12).
11.1.1 PYTHIA interface
The PYTHIA [18] interface is implemented in terms of UPINIT/UPEVNT proce-
dures and 2 → N phase-space generator.
The momenta of external particles are retrieved from PYTHIA user process
event common block (PUP(I,J)). The order of particles in kinematics should
meet the order of generated particles, and the CODE_PYTHIA_VECTORS option
should be used.
The template UPINIT/UPEVNT procedures for ALHEP → PYTHIA junction are
found at ALHEP website. However one should modify them by adjusting the
generated particles sequence, including user cutting rules, using symmetries
for calculation similar processes in single function etc.
The plane 2 → N phase-space generator in alhep_lib.F library file is written
by V. Mossolov. It is desirable to replace the plain phase-space generator by
the adaptive one for multiparticle production process.
The more automation will be implemented in future. Please refer to ALHEP
website for details.
11.2 Mathematica
The output interface to Mathematica Notebook [16] file is basically used to
view the expressions in the convenient form. Implemented for all the symbolic
objects in ALHEP.
SaveNB( (str)fn, (expr or me)val [, (str)comm ][, (int)flags ]);
MarkNB( (str)fn [, (str)comm ][, (int)flags ]);
fn: Mathematica output file name.
val: expression to be stored.
comm: comment text to appear in output file ("" if no comments are required).
flags: file open flags (default is ”0”):
0: append to existing file and do not close it afterward,
FILE_START: delete previous, start new file and add Mathematica header,
FILE_CLOSE: add closing Mathematica block.
The MarkNB function is used to add comments only.
The Mathematica program can open valid files only. The valid x.nb file should
be started (FILE_START) and closed once (FILE_CLOSE). It is convenient to in-
sert MarkNB("x.nb","",FILE_START) and MarkNB("x.nb","",FILE_CLOSE)
calls to start and end of your script file.
No backward interface (Mathematica → ALHEP) is currently available.
11.3 LaTeX
The LaTeX interface in ALHEP is only implemented for Feynman diagrams
drawing. The diagrams are illustrated in terms of AxoDraw [23] package.
DrawDiagrams( (diagrams)d, (str)fn [, (int)flags][, (int)draw]);
MarkTeX( (str)fn, [, (str)comm][, (int)flags]);
d: diagrams set.
fn: output LaTeX file name.
flags: file open flags (0(default): append, FILE_START: truncate old):
comm: comment text to appear in output file ("" if no comments are required).
draw: draw options (default is DD_SMALL):
DD_SMALL: diagrams with small font captions.
DD_LARGE: diagrams with large font captions.
DD_MOMENTA: print particles momenta.
DD_SWAP_TALES: allow arbitrary order for final-state lines (option is included
automatically for 2 → N kinematics).
DD_DONT_SWAP_TALES: deny the DD_SWAP_TALES option. The order of final
lines is same to kinematics declaration.
11.4 ALHEP native save/load operations
The native ALHEP serialization format is XML-structured. It may be edited
outside the ALHEP for some debug purposes.
Save( (str)fn, (diagrams or expr)val);
object Load( (str)fn );
fn: output XML file name.
val: diagrams set or symbolic expression to be stored.
12 Common algebra utilities
The following set of common symbolic operations is available:
• Expand(): expands all the brackets and arranges result. Works slowly with
huge expressions.
• Arrange(): arranges expression (makes alphabetic order in commutative se-
quences). The most of ALHEP functions performs arranging automatically
and there is no need to call Arrange() directly.
• Minimize(): reduces the number of ”sum-multiply” operations in expres-
sion. Should be used for simplification of big sums before numerical calcu-
lations.
• Factor(): factorize expression.
• KinArrange(): arranges expression using kinematic relations,
• KinSimplify(): simplify expression using the kinematic relations. Works
very slowly with large expressions. Mostly useful for 2 → 2 process.
expr Expand( (expr)src );
expr Arrange( (expr)src );
expr Minimize( (expr)src [, (int)flags]);
expr Factor( (expr)src [, (int)flags]);
expr KinArrange( (expr)src [, (int)flags]);
expr KinSimplify( (expr) src [, (int)flags]);
src: source expression (remains unaffected).
Minimize() function flags (default is MIN_DEFAULT):
MIN_DEFAULT = MIN_FUNCTIONS|MIN_DENS|MIN_NUMERATORS,
MIN_FUNCTIONS: factorize functions,
MIN_DENS: factorize denominators,
MIN_NUMERATORS: factorize numerators,
MIN_NUMBERS: factorize numbers,
MIN_ALL_DENOMINATORS: factorize all denominators 1/(a + b) + 1/x →
(x+ a+ b)/(x ∗ (a+ b))),
MIN_ALL_SINGLE_DENS: factorize single denominators (but not sums, prod-
ucts etc.),
MIN_VERIFY: verify result (self-check: expand result back and compare to
source).
Factor() function flags (default is 0):
FACT_NO_NUMBERS: do not factorize numbers,
FACT_NO_DENS: do not factorize fraction denominators,
FACT_ALL_DENS: factorize all denominators,
FACT_ALL_SINGLE_DENS: factorize all single denominators (but not sums,
products etc.),
FACT_VERIFY: verify result (self-check: expand result back and compare to
source).
KinArrange() function flags (default is 0):
KA_MASS_EXACT: do not truncate masses (neglecting SetMassCalcOrder()
settings),
KA_MANDELSTAMS: involve Mandelstam variables (for 2 → 2 kinematics),
KA_NO_EXPAND: do not expand source.
KinSimplify() function flags (default is KS_FACTORIZE_DEFAULT):
KS_FACTORIZE_DEFAULT: factorize functions and denominators (the first
two flags below),
KS_FACTORIZE_FUNCTIONS: factorize functions before simplification,
KS_FACTORIZE_DENS: factorize denominators (including partial factoriza-
tion) before simplification,
KS_FACTORIZE_ALL_DENS: factorize all denominators,
KS_MASS_EXACT: do not truncate masses (neglecting SetMassCalcOrder()
settings),
KS_MANDELSTAMS: involve Mandelstam variables (for 2 → 2 kinematics),
KS_NO_EXPAND: do not expand source (if no simplification are found).
13 Debugging tools
ALHEP allows user to control the most of internal calculation flow. The debug
info is stored to debug.nb file (critical messages will also appear in console).
The debug.ini file contains numerical criteria for messages to be logged, the
debug levels for different internal classes. Warning: raising the debug.ini val-
ues leads to sufficient performance drop and enormous debug.nb file growth.
If one feels some problems with ALHEP usage, please contact author for assis-
tance (attaching your script file). Do not waste the time for manual debugging
using debug.ini.
There are some specific commands to view internal data. For example, the
whole list of tensor virtual integrals reduction results is kept in internal data
storage and can be dumped to Mathematica file for viewing:
• str ViewParticleData((int)PID) returns the brief information on parti-
cle settings. Puts the information string to console and returns it also as a
string variable.
PID: particle ID in kinematics.
• ViewFeynmanRules((str)nb_file, (int)flags) stores Feynman rules of
current physics to Mathematica file.
• ViewTensorVITable((str)nb_file, (int)flags)) stores tensor integrals
reduction table to Mathematica file. The VI reduction table is filled during
the Evaluate() operation.
• ViewScalarVICache((str)nb_file, (int)flags)) stores scalar loop in-
tegrals values cache. The scalar VI cache is filled during the CalcScalarVI()
invocation.
nb file: output Mathematica file name,
flags: access flags to Mathematica file: 0,FILE_START and/or FILE_CLOSE.
Call without parameters turns output to debug.nb file.
14 System commands
The two system commands are useful:
• Halt(): stop further script processing. May be used to test the first part of
script and save (Save()) internal result. It first part finished successfully, it
may be commented (\*...*\) and followed by loading procedure (Load()).
Then script execution is restarted.
• Timer(): view time elapsed since the last Timer() call (from program start
for first call).
15 Examples
15.1 Amplitude for qq̄ → W+W−γ
Let’s consider the example from sec. 2.1 in details. We also extend it for anoma-
lous quartic gauge boson interactions [28]. And we don’t use the Bondarev
method for traces calculation this time. Please refer to sec. 2.1 for ALHEP
installation notes.
We start test.al script from output files creation:
nbfile = "uuWWA_MPXXX.nb"; // Mathenatica file name
MarkNB(nbfile, "", FILE_START); // Create file
texfile = "res.tex"; // LaTeX file name
MarkTeX(texfile, "", FILE_START); // Create file
Then 2 → 3 process kinematics and physics are declared:
u(p1, e1) ū(p2, e2) → γ(f0, g0) W+(f1, g1) W−(f2, g2). (5)
SetKinematics(2, 3 // 2->3 process
,QUARK_U,"p\_1","e\_1" // u
,-QUARK_U,"p\_2","e\_2" // u-bar
,WBOZON,"f\_1","g\_1" // W{+}
,PHOTON,"f\_0","g\_0" // photon
,-WBOZON, "f\_2", "g\_2" ); // W{-}
SetDiagramPhysics(PHYS_SM_Q1GEN|PHYS_4BOSONS_ANOMALOUS);
We declare physics with u- and d-quarks only. The amplitude will be sum-
marized for all the possible internal quarks numerically. It requires the simple
replacing of quark mixing matrix in resulting Fortran code: U2ud → U2ud+U2us+
U2ub).
Next we declare the u- and d-quarks massless:
SetMassCalcOrder(QUARK_U, 0); // consider massless
SetMassCalcOrder(QUARK_D, 0); // consider massless
Set polarizations to ”-+UUU”:
SetFermionHelicity(1, -1); // u
SetFermionHelicity(2, 1); // u-bar
Create diagrams set and store it to LaTeX file:
diags = ComposeDiagrams(3); //e^3 order
DrawDiagrams(diags, texfile);
Next we include the following lines:
Save("diags.xml",diags); //save to XML file
//Halt(); //stop execution
//diags = Load("diags.xml"); //load from XML file
We can save diagrams, stop the program now and view diagrams generated.
To stop ALHEP session the Halt() line should be uncommented. Then we
modify our script as follows:
Fig. 1. The diagrams generated for uū → W−W+γ process (see res.tex file). The
anomalous quartic gauge boson interaction affects the first two diagrams.
/* diags = ComposeDiagrams(3); // commented
... // commented
Halt(); */ // commented
diags = Load("diags.xml"); // uncommented
If we run the script again, it will skip the diagrams generation step and load
diagrams from XML file.
Matrix element retrieval:
me = RetrieveME(diags); //get matrix element
SaveNB(nbfile, me, "Matrix element"); //view
Calculate helicity amplitude, arrange result and minimize the +/× operations
number:
ampl = CalcAmplitude(me);
SaveNB(nbfile, ampl, "Amplitude after CalcAmplitude()");
ampl = KinArrange(ampl);
SaveNB(nbfile, ampl, "Amplitude after KinArrange()");
ampl = Minimize(ampl);
SaveNB(nbfile, ampl, "Amplitude after Minimize()");
The another breakpoint can be inserted here. The result for amplitude is saved,
the Halt and Load commands are commented for further use:
Save("ampl.xml", ampl); // save amplitude
//MarkNB(nbfile, FILE_CLOSE); Halt(); // close NB and exit
//ampl = Load("ampl.xml"); // load amplitude
This breakpoint allows to repeat the next Fortran creation step without re-
calculating of matrix element.
Let’s average over final state polarizations in further numerical procedure. Set
final particles unpolarized:
SetPolarized(-1, 0); // set unpolarized
SetPolarized(-2, 0); // set unpolarized
SetPolarized(-3, 0); // set unpolarized
The Fortran output for differential cross section:
SetParameterS(PAR_STR_FORTRAN_TEMP, "TMP1");
f = NewFortranFile("uuWWA.F", CODE_F77); //f77 file
CreateFortranProc(f, "uuWWA", ampl,
CODE_IS_AMPLITUDE| //square amplitude
CODE_CHECK_DENOMINATORS| //check 1/0 limits
CODE_COMPLEX16| //complex values
CODE_POWER_PARAMS| //F(M^2) instead of F(M)
CODE_PYTHIA_VECTORS); //use PYTHIA PUP(I,J) vectors
The SetParameterS call sets the unique notation for internal variables and
functions. Please do not make it too long. The complex-type code is required
for proper amplitude calculation.
Close Mathematica output file at the end of script:
MarkNB(nbfile, FILE_CLOSE);
The execution of this script takes less than 2 minutes at 1.8GHz P4 processor.
We will not discuss the structure of generated uuWWA.F file in details. But
some remarks should be done:
Line 5: The main function call. The following parameters are declared (order
is changed here): All the parameters (except the kQOrig) are of COMPLEX*16
type. Ones the CODE_COMPLEX16 option is set, all the real objects are treated
as complex.
kQOrig (INTEGER): The ID of u(first)-quark in PYTHIA PUP(I,J) array.
Possible values: 1 or 2.
PAR a 0, PAR a c, PAR a n, PAR ah c, PAR ahat n: Anomalous quartic
gauge boson interaction constants a0, ac, an, âc, ân [28].
PAR CapitalLambda: Scale factor Λ for anomalous interaction [28].
PAR VudP2: Quark mixing matrix element squared |Uud|2. The U2ud +U2us +
U2ub value may be passed to summarize the whole diagrams (neglecting
quarks masses). The numbers for mixing matrix elements may be obtained
using QMIX_VAL(ID1,ID2), QMIX_SQR_SUM(ID) and QMIX_PROD_SUM func-
tions of alhep_lib.F library.
Line 29: Internal COMMON-block with PAR(XX) array. All the scalar cou-
plings and other compound objects are precalculated and stored in PAR(XX).
Lines 49-55: External momenta initialization from PYTHIA PUP(I,J) ar-
ray. The order of external vectors is expected as follows:
PUP(I,1),PUP(I,2): initial particles. If kQOrig=2 the order is backward:
PUP(I,2), PUP(I,1).
PUP(I,3..N): final particles in the same order as in SetKinematics()
call. One should modify this section (or SetKinematics() parameters)
to make the proper particles order.
DO 10 I=1,4
p_1(I) = DCMPLX(PUP(I,kQ1Orig),0D0) // kQ1Orig = kQOrig
f_0(I) = DCMPLX(PUP(I,4),0D0)
p_2(I) = DCMPLX(PUP(I,kQ2Orig),0D0) // kQ2Orig = 3-kQOrig
f_1(I) = DCMPLX(PUP(I,3),0D0)
f_2(I) = DCMPLX(PUP(I,5),0D0)
10 CONTINUE
Lines 69-231: Polarization averaging and basis vector generation cycle. For
any momenta set the PAR(XXX) array is filled. Then denominator checks
and amplitude averaging are performed.
Lines 244-252, 438-444, ... Interaction constants and particle masses defi-
nitions in sub-procedures. The constants can be declared as main function
parameters using CODE_NO_CONSTANTS option in CreateFortranProc func-
tion.
For the complete pp̄ → W+W−γ analysis the following steps are required:
• The PYTHIA client program should be written. The template files are avail-
able at ALHEP website.
• The another helicity configuration +-UUU should be calculated separately.
• The another channels qiq̄j → W+W−γ (i 6= j) should be calculated and
included into generator.
15.2 Z-boxes for e−e+ → µ−µ+
Let’s calculate some box diagrams now. Consider the following process:
e−(p1, e1) e
+(p2, e2) → µ−(f1, g1) µ+(f2, g2). (6)
As in previous example, we start command script from files initialization:
nbfile = "Zbox.nb"; // Mathenatica file name
MarkNB(nbfile, "", FILE_START); // Create file
texfile = "res.tex"; // LaTeX file name
MarkTeX(texfile, "", FILE_START); // Create file
The 2 → 2 kinematics declaration:
SetKinematics(2, 2, // 2->2 process
ELECTRON, "p\_1", "e\_1" , // e^{-}
-ELECTRON, "p\_2", "e\_2", // e^{+}
MUON, "f\_1", "g\_1", // mu^{-}
-MUON, "f\_2", "g\_2"); // mu^{+}
The Standard model physics in Feynman gauge (we omit quarks for faster
diagrams generation):
SetDiagramPhysics(PHYS_ELW|PHYS_NOQUARKS);
Set leptons massless and declare the N -dimensional space:
SetMassCalcOrder(ELECTRON, 0); // massless electrons
SetMassCalcOrder(MUON, 0); // massless muons
SetNDimensionSpace(1);
Use Mandelstam variables throughout the calculation:
SetParameter(PAR_MANDELSTAMS,1);
Consider unpolarized particles:
SetPolarized(1, 0); // unpolarized e^{-}
SetPolarized(2, 0); // unpolarized e^{+}
SetPolarized(-1, 0); // unpolarized mu^{-}
SetPolarized(-2, 0); // unpolarized mu^{+}
Compose born and one-loop diagrams:
diags_born = ComposeDiagrams(2); //e^2 order
DrawDiagrams(diags_born, texfile);
diags_loop = ComposeDiagrams(4); //e^4 order
Save("diags_loop.xml",diags_loop);
//diags_loop=Load("diags_loop.xml");
DrawDiagrams(diags_loop, texfile);
The 220 loop diagrams are created, saved to internal format and TeX file.
The ComposeDiagrams(4) procedure takes 5-10 minutes here. It is convenient
to comment ComposeDiagrams(4)-Save() lines at second run and Load loop
diagrams from disk.
Fig. 2. Born level diagrams for e−e+ → µ−µ+.
Fig. 3. Part of 220 loop diagrams stored to res.tex file
Let’s select the double Z-exchange box graphs from the whole set (194 and
195 diagrams at fig. 3):
diag_box = SelectDiagrams(diags_loop,194,195);
Next we couple the loop and born matrix elements:
me_born = RetrieveME(diags_born);
me_box = RetrieveME(diag_box);
me_sqr = SquareME(me_box, me_born);
The simplification procedures are not included into SquareME implementation.
It may take much time to arrange items in huge expression. Therefore all the
simplification procedures are optional and should be called manually:
me_sqr = KinArrange(me_sqr);
me_sqr = KinSimplify(me_sqr);
SaveNB(nbfile, me_sqr, "squared & simplified");
The reduction of tensor virtual integrals follows:
me_sqr = Evaluate(me_sqr);
me_sqr = KinArrange(me_sqr);
me_sqr = KinSimplify(me_sqr);
SaveNB(nbfile, me_sqr, "VI evaluated");
Next we convert scalar integrals to invariant-dependent form and replace with
tabulated values:
me_sqr = ConvertInvariantVI(me_sqr);
me_sqr = CalcScalarVI(me_sqr); // use pre-calculated values
me_sqr = KinArrange(me_sqr);
SaveNB(nbfile, me_sqr, "VI scalars ");
Turn to 4-dimensional space, drop out (n− 4)i items and final simplification:
me_sqr = SingularArrange(me_sqr);
SetNDimensionSpace(0);
me_sqr = KinArrange(me_sqr);
me_sqr = KinSimplify(me_sqr);
Save result and create Fortran code with LoopTools [19] interface:
Save("ZBox.xml",me_sqr); // save result
//me_sqr = Load("ZBox.xml"); // reload result
SaveNB(nbfile, me_sqr, "Z boxes result"); // view result
f = NewFortranFile("ZBOX.F", CODE_F77);
CreateFortranProc(f, "ZBOX", me_sqr,
CODE_POWER_PARAMS|CODE_LOOPTOOLS);
View tensor integrals reduction table and close Mathematica output file:
ViewTensorVITable(nbfile);
MarkNB(nbfile, FILE_CLOSE);
The script runs about 15 minutes on 1.8GHz P4 processor. The half of this
time takes the ComposeDiagrams(4) procedure.
�������
jjjjjjjj
jjjjjjjj
2,t,me
2,s,me
2,mZ,m�,me,mZ<
jjjjjjjj
jjjjjjjj
jjjjjjjj
jjjjjjj
����������������������HmZ2 − sL JmZ
3N s − JmZ2 t2 + mZ4 t +
3N e2
zzzzzzz gm
����������������������HmZ2 − sL s t
zzzzzzzz
zzzzzzzz
jjjjjjjj
����������������������HmZ2 − sL JmZ
3N s − JmZ2 t2 + mZ4 t + 1����
3N e2
zzzzzzzz
zzzzzzzz
2,s,m
,mZ,mZ<
jjjjjjjj
jjjjjjjj
J 1����
����������������������HmZ2 − sL s − J
uN e2
zzzzzzzz
jjjjjjj
jjjjjjj
jjjjjjjJ
����������������������HmZ2 − sL s − J
uN e2
zzzzzzz gm
HmZ2 + 1����2 s + uL
������������������������������������HmZ2 − sL s
zzzzzzz gm
2 JmZ2 + 1����
s + uN e2
zzzzzzz gm
zzzzzzzz
2,s,me
2,me,mZ,mZ<
jjjjjjjj
jjjjjjjj
J 1����
����������������������HmZ2 − sL s − J
uN e2
zzzzzzzz
jjjjjjj
jjjjjjj
jjjjjjjJ
����������������������HmZ2 − sL s − J
uN e2
zzzzzzz gm
HmZ2 + 1����2 s + uL
������������������������������������HmZ2 − sL s
zzzzzzz gm
2 JmZ2 + 1����
s + uN e2
zzzzzzz gm
zzzzzzzz
s + D0
2,u,me
2,s,me
2,mZ,m�,me,mZ<
jjjjjjjj
jjjjjjj
jjjjjjJ2 mZ
t + 2 mZ
����������������������HmZ2 − sL s − J2 mZ
t + 2 mZ
3N e2
zzzzzz
jjjjjjj
����������������������HmZ2 − sL s − e
zzzzzzz gm
zzzzzzz gm
jjjjjjjj
����������������������HmZ2 − sL s − e
zzzzzzzz
zzzzzzzz
2,t,me
2,mZ,m
HmZ2 + 1����2 t − 1����2 uL
�������������������������������������������HmZ2 − sL s − JmZ
uN e2y
zzz g
HmZ2 + 1����2 t − 1����2 uL
�������������������������������������������HmZ2 − sL s − JmZ
uN e2y
zzz g
zzz t −
8s,mZ,mZ<
jjjjjjjj
jjjjjjjj
jjjjjjjj
jjjjjjj
����������������������HmZ2 − sL s − e
zzzzzzz gm
����������������������HmZ2 − sL s t
zzzzzzzz
zzzzzzzz
jjjjjjjj
����������������������HmZ2 − sL s − e
zzzzzzzz
zzzzzzzz
2,u,me
2,mZ,m�,me< i
jjjjjjH2 mZ
+ s + 2 uL g
����������������������HmZ2 − sL s − H2 mZ
+ s + 2 uL e2
zzzzzz g
LogA 1�������
jjjjjjjj
jjjjjjj
����������������������HmZ2 − sL s − e
zzzzzzz gm
jjjjjjjj
����������������������HmZ2 − sL s − e
zzzzzzzz
zzzzzzzz
LogA 1�������
jjjjjjg
����������������������HmZ2 − sL s − e
zzzzzz g
jjjjjjjj
H2 + χL g
jjjjjjjj
H2 + χL
����������������������HmZ2 − sL s t −
jjjjjjjJ1 +
����������������������HmZ2 − sL s − J1 +
χN e2
zzzzzzz gm
zzzzzzzz
zzzzzzzz
jjjjjjjj
J1 + 1����
����������������������HmZ2 − sL s + e
2 J−J1 + 1����
zzzzzzzz
zzzzzzzz
zzzzzzzz
Fig. 4. The result expression in Zbox.nb file. This result contains UV-regulator
term χ, that should cancel one in B0(s,MZ ,MZ) integral. It can be checked using
GetUVTerm() function.
Code remarks:
Line 5: The main function call. The 3 parameters are usual Mandelstam
variables (s, t, u, type is complex). The current ALHEP version does not
care about interdependent parameters in Fortran output. And all the three
Mandelstam variables may occur in parameters list. The future versions will
be saved from this trouble.
Line 14,50: Include LoopTools header file ("looptools.h"). See LoopTools
manual [19] for details.
Line 20,21: Retrieve LoopTools values for UV-regulator getdelta() and
DR-mass squared getmudim().
The complete code of examples including scripts, batches and output files are
available at ALHEP website.
16 Conclusions
The new program for symbolic computations in high-energy physics is pre-
sented. In spite of several restrictions remained in current version, it can be
useful for computation of observables in particle collision experiments. It con-
cerns both multiparticle production amplitudes and loop diagrams analysis.
The nearest projects are:
• Bondarev functions method improvement,
• Complete renormalization scheme for SM,
• Complete covariant analysis of the one-loop radiative corrections including
the hard bremsstrahlung scattering contribution.
• Arbitrary Lagrangian assignment.
Refer ALHEP project websites for program updates:
http://www.hep.by/alhep ,
http://cern.ch/~makarenko/alhep .
References
[1] CompHEP, http://theory.sinp.msu.ru/comphep,
E. Boos et al., Nucl.Instrum.Meth. A534 (2004) 250, hep-ph/0403113.
[2] SANC, http://sanc.jinr.ru/ (or http://pcphsanc.cern.ch/),
A. Andonov et al., hep-ph/0411186, to appear in Comp.Phys.Comm.; D.
Bardin, P. Christova, L. Kalinovskaya, Nucl.Phys.Proc.Suppl. B116 (2003) 48.
[3] GRACE, http://minami-home.kek.jp/,
G. Belanger et al., LAPTH-982/03, KEK-CP-138, hep-ph/0308080 (One-
loop review); J. Fujimoto et al., Comput. Phys. Commun. 153 (2003) 106,
hep-ph/0208036 (SUSY review);
http://theory.sinp.msu.ru/comphep
http://arxiv.org/abs/hep-ph/0403113
http://sanc.jinr.ru/
http://arxiv.org/abs/hep-ph/0411186
http://minami-home.kek.jp/
http://arxiv.org/abs/hep-ph/0308080
http://arxiv.org/abs/hep-ph/0208036
[4] MadGraph, http://madgraph.hep.uiuc.edu/,
T. Stelzer, W. F. Long, Comput. Phys. Commun. 81 (1994) 357,
hep-ph/9401258.
[5] O’Mega, http://theorie.physik.uni-wuerzburg.de/~ohl/omega/
M. Moretti, T. Ohl, J. Reuter, IKDA 2001/06, LC-TOOL-2001-040,
hep-ph/0102195.
[6] FormCalc, http://www.feynarts.de/formcalc,
T. Hahn, M. Perez-Victoria, Comput.Phys.Commun. 118 (1999) 153,
hep-ph/9807565.
[7] FeynCalc, http://www.feyncalc.org/
R. Mertig, M. Bohm, A. Denner, Comput.Phys.Commun. 64 (1991) 345.
[8] Amegic, F. Krauss, R. Kuhn, G. Soff, JHEP 0202 (2002) 044, hep-ph/0109036.
[9] AlpGen, http://mlm.home.cern.ch/mlm/alpgen/,
M. Mangano et al., CERN-TH-2002-129, JHEP 0307 (2003) 001,
hep-ph/0206293.
[10] HELAC-PHEGAS, C. Papadopoulos, Comput.Phys.Commun. 137 (2001) 247,
hep-ph/0007335; A. Kanaki, C. Papadopoulos, hep-ph/0012004.
[11] xloops, http://wwwthep.physik.uni-mainz.de/~xloops/,
L. Brucher, J. Franzkowski, D. Kreimer, Comput. Phys. Commun. 115 (1998)
[12] aITALC, http://www-zeuthen.desy.de/theory/aitalc/,
A. Lorca, T. Riemann, DESY 04-110, SFB/CPP-04-22, hep-ph/0407149;
J. Fleischer, A. Lorca, T. Riemann, DESY 04-161, SFB/CPP-04-38,
hep-ph/0409034.
[13] MINCER,
http://www.nikhef.nl/~form/maindir/packages/mincer,
S.A. Larin, F.V. Tkachov, J.A.M. Vermaseren, NIKHEF-H/91-18;
S.G. Gorishnii, S.A. Larin, L.R. Surguladze, F.V. Tkachov, Comput. Phys.
Commun. 55 (1989) 381.
[14] DIANA, http://www.physik.uni-bielefeld.de/~tentukov/diana.html,
M. Tentyukov, J. Fleischer, Comput.Phys.Commun. 160 (2004) 167,
hep-ph/0311111; M. Tentyukov, J. Fleischer, Comput.Phys.Commun. 132
(2000) 124, hep-ph/9904258.
[15] REDUCE by A. Hearn, http://www.reduce-algebra.com/.
[16] Mathematica by S.Wolfram,
http://www.wolfram.com/products/mathematica/.
[17] FORM by J. Vermaseren, http://www.nikhef.nl/~form/.
[18] PYTHIA 6.4, http://projects.hepforge.org/pythia6/,
T.Sjostrand, S.Mrenna, P.Skands, JHEP 0605 (2006) 026, LU-TP-06-13,
hep-ph/0603175.
http://madgraph.hep.uiuc.edu/
http://arxiv.org/abs/hep-ph/9401258
http://arxiv.org/abs/hep-ph/0102195
http://arxiv.org/abs/hep-ph/9807565
http://www.feyncalc.org/
http://arxiv.org/abs/hep-ph/0109036
http://mlm.home.cern.ch/mlm/alpgen/
http://arxiv.org/abs/hep-ph/0206293
http://arxiv.org/abs/hep-ph/0007335
http://arxiv.org/abs/hep-ph/0012004
http://www-zeuthen.desy.de/theory/aitalc/
http://arxiv.org/abs/hep-ph/0407149
http://arxiv.org/abs/hep-ph/0409034
http://arxiv.org/abs/hep-ph/0311111
http://arxiv.org/abs/hep-ph/9904258
http://www.reduce-algebra.com/
http://arxiv.org/abs/hep-ph/0603175
[19] LoopTools, http://www.feynarts.de/looptools/,
T. Hahn, M. Perez-Victoria, Comput.Phys.Commun. 118
(1999) 153, hep-ph/9807565; T. Hahn, Nucl.Phys.Proc.Suppl. 89 (2000) 231,
hep-ph/0005029.
[20] FF, http://www.xs4all.nl/~gjvo/FF.html,
G.J. van Oldenborgh, J.A.M. Vermaseren, Z.Phys. C46 (1990) 425, NIKHEF-
H/89-17.
[21] FeynArts, http://www.feynarts.de/,
T. Hahn, Comput.Phys.Commun. 140 (2001) 418, hep-ph/0012260; T. Hahn, C.
Schappacher, Comput.Phys.Commun. 143 (2002) 54, hep-ph/0105349 (MSSM).
[22] WHIZARD, http://www-ttp.physik.uni-karlsruhe.de/whizard/,
W. Kilian, LC-TOOL-2001-039.
[23] AxoDraw LaTeX style package,
J.A.M. Vermaseren, Comput.Phys.Commun. 83 (1994) 45.
[24] A.L. Bondarev, Talk at 9th Annual RDMS CMS Collaboration Conference,
Minsk, 2004, hep-ph/0511324.
[25] A.L. Bondarev, Nucl.Phys. B733 (2006) 48, hep-ph/0504223.
[26] P. De Causmaecker, R. Gastmans, W. Troost, Tai Tsun Wu, Nucl. Phys. B206
(1982) 53.
[27] A. Denner, Fortsch.Phys. 41 (1993) 307.
[28] A.Denner, S.Dittmaier, M.Roth,D.Wackeroth, Eur.Phys.J.C 20 (2001) 201.
http://www.feynarts.de/looptools/
http://arxiv.org/abs/hep-ph/9807565
http://arxiv.org/abs/hep-ph/0005029
http://www.feynarts.de/
http://arxiv.org/abs/hep-ph/0012260
http://arxiv.org/abs/hep-ph/0105349
http://www-ttp.physik.uni-karlsruhe.de/whizard/
http://arxiv.org/abs/hep-ph/0511324
http://arxiv.org/abs/hep-ph/0504223
	Introduction
	ALHEP Review. Program Structure 
	Getting Started 
	Calculation scheme 
	ALHEP script language 
	Initialization section 
	Particles 
	Kinematic selection 
	Particle masses 
	Polarization data 
	Physics selection 
	Bondarev functions 
	Diagrams generation 
	Helicity amplitudes 
	Matrix Element squaring 
	Virtual integrals operations 
	Regularization 
	ALHEP interfaces 
	Fortran numerical code 
	Mathematica
	LaTeX 
	ALHEP native save/load operations 
	Common algebra utilities 
	Debugging tools 
	System commands 
	Examples 
	Amplitude for q "7016q W+ W-  
	Z-boxes for e- e+ - + 
	Conclusions
	References
