What is xloops?

xloops was written to implement the methods developed in Mainz for calculating one- and two-loop Feynman diagrams. Whereas many other packages for loop calculations are often ad-hoc implementations computing numeric results for very specific processes or special cases (but quite efficiently so), the idea of xloops is to have a program that works for arbitrary massive diagrams to any tensor rank, and that can also return analytic results if possible. Eventually, xloops will be expanded into a package for the complete automatic calculation of processes in the Standard Model.

The roots of xloops go back to a Maple package for one-loop integrals released in 1993. Due to inherent problems with Maple as an implementation environment, the project was switched to C++ in 1999. One result of this was the creation of GiNaC, a C++ library for symbolic computation, that is now a project of its own.

Currently, xloops consists of a library of C++ and GiNaC functions for calculating arbitrary one-loop one-, two-, and three-point, and two-loop self-energy integrals of any tensor rank. All one-loop integrals and the divergent part of two-loop integrals are returned in symbolic form as an expansion in the dimensional regularization parameter ε. The finite part of two-loop integrals is transformed into a two-dimensional function that is numerically integrable.

Examples

  1. A one-loop calculation
    Here is a complete example program that calculates the divergent and finite parts of the general rank-1 one-loop two-point integral: 
    #include <iostream>
#include <ginac/ginac.h>
#include <xloops/xloops.h>

using std::cout; using std::endl;
using namespace GiNaC;
using namespace xloops;

int main(void)
{
    // Declare symbols for momentum, masses, and parameters
    symbol q("q"), m1("m1"), m2("m2"), eps("eps"), rho("rho");

    // One-loop integral function
    ex C = OneLoop2Pt(1, 0, q, m1, m2, 1, 1, rho, eps);

    // Perform epsilon-expansion
    C = C.series(eps == 0, 1);

    // Extract the divergent and finite parts
    ex div = C.coeff(eps, -1).subs(rho == 0).normal();
    ex fin = C.coeff(eps, 0).subs(rho == 0);

    // Output the result
    cout << "Order eps^-1 is " << endl << div << endl;
    cout << "Order eps^0 is "  << endl << fin << endl;
    return 0;
}
  2. A two-loop calculation
    The following program first computes the divergent part of the scalar two-loop sunrise topology. Then it evaluates the finite part of the integral for a specific choice of the masses and the external momentum. 
    #include <iostream>
#include <ginac/ginac.h>
#include <xloops/xloops.h>

using std::cout; using std::endl;
using namespace GiNaC;
using namespace xloops;

int main(void)
{
    // Declare symbols for momentum, masses, and parameters
    symbol q("q");
    symbol m1("m1"), m2("m2"), m3("m3");
    symbol rho("rho"), eps("eps");
    symbol x("x"), y("y");

    // Two-loop two-point integral function, sunrise topology
    ex C = TwoLoop2Pt4(0, 0, 0, 0, 0, q, lst(m1, m2, m3), lst(1, 1, 1), rho, eps, lst(x, y));

    // The result is a two-element list containing the analytic and numeric
    // parts; extract them
    ex num = C.op(0);
    ex ana = C.op(1);

    // Perform epsilon-expansion of the analytic part
    ana = ana.series(eps == 0, 1);

    // Extract the divergent (eps^-2 and eps^-1) and finite parts
    ex eps2 = collect_common_factors(ana.coeff(eps, -2).subs(rho == 0).normal());
    ex eps1 = collect_common_factors(ana.coeff(eps, -1).subs(rho == 0).normal());
    ex eps0 = ana.coeff(eps, 0);

    // Output the result for the divergent part
    cout << "eps^-2 = " << eps2 << endl;
    cout << "eps^-1 = " << eps1 << endl;

    // Substitute values for masses and momentum in the finite parts
    lst s;
    s = rho == 1E-15,
        m1 == 2,
        m2 == 4,
        m3 == 6,
        m4 == 8,
        m5 == 10,
        q == 20;
    eps0 = eps0.subs(s).evalf();
    fin = fin.subs(s);

    // Perform numerical integration using Vegas
    fin = VNumInt(fin, 10000, 10, 100000, 10, lst(x, y));

    // Output the result for the finite part, including the finite
    // analytic term
    cout << "finite = " << fin.op(0) + fin.op(2)*I + eps0 << endl;
    cout << "error real part = " << fin.op(1) << endl;
    cout << "error imag part = " << fin.op(3) << endl;

    return 0;
}