This is the program used to calculate the complete next-to-leading cross section for inclusive hadron production in pA collisions, described in the paper

Anna M. Stasto, Bo-Wen Xiao, David Zaslavsky
"Towards the Test of Saturation Physics Beyond Leading Logarithm"
arXiv:1307.4057 [hep-ph]

Please cite this paper if you use the results of the code in a publication.

Quick start

Downloading the code

If you have git installed

The quickest way to obtain the code is to use git. Run the following commands:

git clone https://github.com/diazona/SOLO.git
cd SOLO
git submodule init
git submodule update

If you don't have git installed

If you don't have git (and don't want to download it), you can download the source code manually.

Download the latest designated release as a GZIP'd TAR file or a ZIP file. It contains the following directory structure:
- SOLO-0.1/
  - cmake/ (contains files)
  - interp2d/ (empty)
  - quasimontecarlo/ (empty)
  - src/ (contains files)
  - README
  - CMakeLists.txt
  - .gitmodules
Extract this into a folder of your choice.
Download the latest release of interp2d as a GZIP'd TAR file or a ZIP file. It contains the following directory structure:
- interp2d-0.1/
  - CMakeLists.txt
  - COPYING
  - Doxyfile.in
  - README.md
  - bicubic.c
  - bilinear.c
  - interp2d.c
  - interp2d.h
  - interp2d.kdev4
  - interp2d_spline.c
  - interp2d_spline.h
  - test.c
Move the files from interp2d-0.1/ into SOLO-0.1/interp2d/.
Download the latest release of quasimontecarlo as a GZIP'd TAR file or a ZIP file. It contains the following directory structure:
- quasimontecarlo-0.1/
  - CMakeLists.txt
  - README.md
  - cmake/ (contains files)
  - quasimontecarlo.c
  - quasimontecarlo.h
  - quasimontecarlo.kdev4
  - test.c
Move the files from quasimontecarlo-0.1/ into SOLO-0.1/quasimontecarlo/.

Downloading external files

At this point you will need to manually download the MSTW PDF interface from http://mstwpdf.hepforge.org/code/code.html. Extract the files mstwpdf.cc and mstwpdf.h from the tarball and place them in SOLO's src/ directory. (We are not authorized to distribute the MSTW PDF interface as part of SOLO, which is why this has to be done manually.)

Compiling

You will need to have the following dependencies already installed:

A C++ compiler and runtime environment (GCC and glibc are recommended)
GSL, the GNU Scientific Library
OpenSSL

If you have CMake installed

CMake is a cross-platform build system that fills roughly the same role as configure from GNU autotools. If you have it installed, it will check that the dependencies are satisfied, and show an error message if they are not.

To compile with CMake: from the directory SOLO/ or SOLO-0.1/ (the parent directory of src/), run these commands:

mkdir build
cd build
cmake .. && make

At the end of this you should have a build/ directory containing the executable files.

If you don't have CMake installed

If you don't have CMake installed (and don't want to download it), you can compile the source code manually, but it is up to you to make sure all the dependencies are properly installed. In the directory src/, make an empty file named git_revision.h. Then, from the parent directory of src/, run this command:


c++ src/oneloopcalc.cpp src/integrator.cpp src/hardfactors_position.cpp \
src/hardfactors_radial.cpp src/hardfactors_momentum.cpp src/dss_pinlo.cpp \
src/mstwpdf.cc src/cubature.c src/gluondist.cpp src/coupling.cpp \
src/factorizationscale.cpp src/integrationcontext.cpp \
src/integrationtype.cpp src/context.cpp src/utils.cpp interp2d/interp2d.c \
interp2d/bilinear.c interp2d/bicubic.c quasimontecarlo/quasimontecarlo.c \
-I interp2d/ -I quasimontecarlo/ -lgsl -lgslcblas -lm -lssl -lcrypto \
-o oneloopcalc

At least, this is the command for GCC; if you are using a different compiler, adjust it accordingly.

Extra runtime files

In order to run the program, you will need two additional files:

The grid file for the MSTW 2008 PDF at NLO, from the paper

A. D. Martin, W. J. Stirling, R. S. Thorne and G. Watt,
"Parton distributions for the LHC",
Eur. Phys. J. C 63 (2009) 189-285
arXiv:0901.0002 [hep-ph].

The filename is mstw2008nlo.00.dat, and it can be downloaded as part of an archive at the MSTW PDF site.
The data file for the DSS fragmentation functions at NLO, from the paper

Daniel de Florian, Rodolfo Sassot, Marco Stratmann
"Global analysis of fragmentation functions for pions and kaons and their uncertainties"
Phys. Rev. D 75, 114010 (2007)
arXiv:hep-ph/0703242

The filename is PINLO.DAT, and it is reportedly available on request from the authors of the DSS paper. Unfortunately we (authors of SOLO) are not aware of a website where this file is directly downloadable.

Running

The main executable is called oneloopcalc. It should be invoked as

./oneloopcalc <options>

where the <options> can include any number of the following, in any order:

Hard factor group specifications

Hard factor specifications and hard factor group specifications tell the program which terms to calculate out of the complete NLO expression.

A hard factor specification is a string like "p:h02qq" or "m:h16gg". The "p:" at the beginning specifies the position space implementation, "r:" specifies a position space implementation with the angular integral already done, and the "m:" specifies the momentum space implementation. The prefix can be omitted, in which case position space is taken as the default. (Not recommended, as position space is highly inaccurate for some terms.)

The rest of the string gives the name of a hard factor. The set of possible names that can be used with a "p:" prefix are all the return values from the get_name() methods in hardfactors_position.h, and similarly for "r:" with hardfactors_radial.h and "m:" with hardfactors_momentum.h. Here's a list transcribed from those files (but check the source code because this could be out of date):

p:h02qq m:h02qq p:h12qq r:h12qq p:h14qq m:h14qq p:h02gg m:h02gg p:h12gg r:h12gg p:h12qqbar m:h12qqbar p:h16gg m:h16gg p:h112qg r:h112qg p:h122qg r:h122qg p:h14qg m:h14qg p:h112gq r:h112gq p:h122gq r:h122gq p:h14gq m:h14gq

The names are case-insensitive.
A hard factor group specification is made of any number of individual hard factor specifications separated by commas. The program will report only the total of all terms in the group. (A hard factor specification as in the last paragraph can be considered a group of one.)

There are some hard factor groups that have abbreviated names:
- lo is equivalent to m:h02qq,m:h02gg
- nlo is equivalent to r:h12qq,m:h14qq,r:h12gg,m:h12qqbar,m:h16gg,r:h112gq, r:h122gq,m:h14gq,r:h112qg,r:h122qg,m:h14qg
These are defined in oneloopcalc.cpp.

If no hard factor groups are specified on the command line, the program behaves as if lo nlo was given, which prints out the total of the LO terms and the total of the NLO terms, as well as the overall total.

Note that specifying several hard factors as a group (with commas) is considerably more efficient than specifying them individually and having the program print out the results for each one separately.

Configuration file names

Command-line arguments can also include the names of any number of configuration files, which contain parameters for the program in the format

key1=value1
key2=value2

and so on. Keys are case-insensitive. All the possible keys are defined throughout the file context.cpp. Here's a list which should be complete:

A (no default): the mass number of the target nucleus
abserr (default 10^-20): the absolute error at which to stop an integration, for integration strategies which use this termination condition
alphas (default 0.2): if fixed coupling is being used, value for the coupling
beta (default 11-2*Nf/3): if running coupling is being used, coefficient for the LO running coupling
c (no default): the centrality coefficient, in the range 0 (central) to 1 (peripheral)
c0r_optimization (default true): if the factorization scale scheme is c0r, whether to skip calculating terms that should be zero
CF (default 1.5): the color factor
coupling_type (default fixed): fixed or running
cubature_iterations (default 1000000): number of calls to use for cubature integration
exact_kinematics (default false): whether to use exact kinematic expressions
factorization_scale (default fixed): "fixed" to set the factorization scale to a fixed value, or "4pT2" or "CpT2" to set it to a multiple of hadron transverse momentum, or "c0r" to set it to be inversely proportional to transverse position
factorization_scale_coefficient (no default): if factorization_scale is "CpT2", this is the coefficient to multiply by pT2 to get mu2 (setting this to 4 is equivalent to specifying factorization_scale = 4pT2)
ff_filename (default PINLO.DAT): filename to read DSS fragmentation function data from
gammaMV (default 1): if gdist_type is set to MV, the anomalous dimension in the MV gluon distribution
gdist_momentum_filename (no default): file to read the momentum-space gluon distribution from, which should contain three columns of space-separated numbers: Q^2, ln(1/x), and F(x,Q^2)
gdist_position_filename (no default): file to read the position-space gluon distribution from, which should contain three columns of space-separated numbers: r^2, ln(1/x), and S(x,r^2)
gdist_subinterval_limit (default 10000): if gdist_type is set to MV, the number of subdivisions to use when performing the Fourier integral of the position-space gluon distribution
gdist_type (default GBW): the type of the gluon distribution: GBW, MV, fMV, file, GBW+file (which uses the GBW distribution at large x and data from a file at small x), or MV+file (same but with MV)
hadron (no default): the type of hadron to be detected, pi-, pi0, or pi+ (this influences which fragmentation function is used)
inf (default 40): the cutoff used for integration over an infinite region (for most integration types, integration over the plane is replaced by integration over a square of half-side equal to the value of inf
integration_strategy (default VEGAS): the integration type to use for Monte Carlo integrals: MISER, VEGAS (best), or QUASI
lambda (default 0.288): if gdist_type is GBW, MV, GBW+file, or MV+file, the exponent in the definition of the saturation scale
lambdaMV (default 0.241): if gdist_type is MV or MV+file, the scale for 1/r in the MV gluon distribution, in GeV
lambdaQCD (default 0.2428711 = sqrt(0.0588)): if coupling_type is running, the QCD scale in GeV
miser_iterations (default 10000000): if integration_strategy is MISER, the number of iterations to use in MISER integration
mu2 (default 10): if factorization_scale is fixed, the factorization scale in GeV
Nc (default 3): number of colors
Nf (default 3): number of flavors
pdf_filename (default mstw2008nlo.00.dat): filename to read MSTW PDF from
projectile (no default): the type of projectile, deuteron or proton
pseudorandom_generator_seed (default 0): seed for the GSL random number generator (0 tells GSL to use a system-dependent default; see the GSL documentation for details)
pseudorandom_generator_type (default mt19937): algorithm to use for generating random numbers; allowed values are in the GSL documentation
pT (no default): comma-separated list of transverse momenta
quasirandom_generator_type (default halten): algorithm to use for generating quasirandom numbers for QMC integration; allowed values are in the GSL documentation
quasi_iterations (default 1000000): the number of iterations at which to stop quasi Monte Carlo integration
regulator (default 1): if coupling is running, the position of the Landau pole for the regulated LO running coupling
relerr (default 0): the relative error at which to stop an integration, for strategies which use this termination condition
satscale_source (default extract from momentum): for a file gluon distribution, how to extract the saturation scale; allowed values are analytic (Q0²(x0/x)^λ), extract from momentum, which determines the saturation scale by finding the momentum where the gluon distribution equals a fixed fraction of its value at a reference momentum, and extract from position, which finds the radius where the gluon distribution equals a fixed threshold value
satscale_threshold (no default): if satscale_source is extract from momentum or extract from position, this is the fixed threshold value (or fraction of its value at a reference point, in the momentum case) that the gluon distribution should equal at the saturation scale
Sperp (default 1): cross-sectional area of the hadron
sqs (no default): sqrt(s), the collider's CM energy
TR (default 0.5): group coefficient
vegas_incremental_iterations (default 1000000): number of function evaluations to use in each step of the VEGAS Monte Carlo algorithm after the first
vegas_initial_iterations (default 100000): number of function evaluations to use to refine the grid in the first step of the VEGAS algorithm
x0 (default 0.000304): the fit parameter from the definition of the saturation scale
Y (no default): comma-separated list of rapidities (in the center of mass frame) to run the calculation at

The configuration files have to at least set A, c, sqs, and Y, and also pT if no transverse momenta are specified as command line arguments.

Transverse momentum values

Any numbers given as command line arguments are put together into one big list of transverse momentum values to run the calculation at. If a comma-separated list of numbers is given, then it will be split apart and each number added to the one big list. There's no significance to putting certain pT values together and others not. (0.5 0.7 0.8,0.9 and 0.5,0.7 0.8 0.9 are exactly equivalent.) Any pT values specified on the command line will replace pT values specified in configuration files.

Literal options

--separate: Print out the results for each individual hard factor, not just the total for each hard factor groups
--minmax: Track and print out the minimum and maximum values of kinematic variables
--trace-gdist: Print out parameters and values for every call to the gluon distribution. Output goes to the file trace_gdist.output in the working directory. (Expect this file to grow to several hundred megabytes.)
--trace=var1,var2,...: Print out selected variables from the integration context after every single evaluation of the function. The output goes to the file trace.output in the working directory. (Expect this file to grow to several megabytes.) The allowable variables are those in ictx_var_list.inc, or you can use --trace=all or --trace=* to print out all available variables.

Structure

Basically the code runs as follows:

Collect the command line options and settings from configuration files and put everything into a ResultsCalculator
For each combination of pT and Y, and for each hard factor group:
1. Create an Integrator with the current values of pT and Y and the current hard factor group
2. The Integrator calls the GSL Monte Carlo integration routine
3. For each time the MC routine evaluates the function
  1. Update the variables in the IntegrationContext
  2. Go through the list of HardFactor instances in the current group and get a value from each one
  3. Return the total value
4. Store the value and error bound returned from the Monte Carlo
Print out all the results

Files

Source code for the program itself

oneloopcalc.cpp: Main program
log.h: Declares an output stream to write status messages to
gsl_exception.h: Declares an exception to be thrown when GSL reports an error
hardfactors.h: A class that abstractly represents a hard factor (i.e. an expression to be integrated)
hardfactors_momentum.h hardfactors_momentum.cpp: Implementation of the momentum space hard factors (terms)
hardfactors_position.h hardfactors_position.cpp: Implementation of the position space hard factors (terms)
gluondist.h gluondist.cpp: Implementations of the gluon distributions
gluondist_driver.cpp: A program to print out values from the gluon distributions
coupling.h coupling.cpp: Implementations of the fixed and LO running couplings
dss_pinlo.h dss_pinlo.cpp: A C++ interface to the DSS fragmentation functions
factorizationscale.h factorizationscale.cpp: Implementations of the various schemes for the factorization scale
integrationcontext.h integrationcontext.cpp: A class that stores the kinematic variables used in the calculation. The values stored in this get updated every time the function is evaluated.
integrator.h integrator.cpp: A class that stores the parameters for the integral and actually calls the GSL Monte Carlo integration functions
integrationtype.h integrationtype.cpp: Definitions of integration types. An integration type specifies how many dimensions are in the Monte Carlo integral and what the limits are.
utils.h utils.cpp: Some string and list processing functions
ictx_var_list.inc: Variables from the integration context, listed in a separate file as a preprocessor hack of sorts

Source code provided separately

interp2d.h libinterp2d.a: A 2D interpolation library compatible with the GSL. Full source code at https://github.com/diazona/interp2d
quasimontecarlo.h libquasimontecarlo.a: A library for quasi Monte Carlo integration compatible with the GSL. Full source code at https://github.com/diazona/quasimontecarlo
mstwpdf.h mstwpdf.cc: A C++ interface to the MSTW PDFs
cubature.h cubature.c: A library for multidimensional cubature (deterministic integration)

Other files

CMakeFiles.txt: Instructions for the build system, CMake
PINLO.DAT: DSS fragmentation function data
pinlo_extended.dat: DSS fragmentation function data with extrapolation to lower z (the significance of this file will be described in an upcoming long paper about SOLO)
mstw2008nlo.00.dat: MSTW PDF data