Global Usage

Design Philosophy

The global category covers the system of units, constants, numerics and random number handling. It can be considered a place-holder for “general purpose” classes used by all categories defined in Geant4. No back-dependencies to other Geant4 categories affect the “global” domain. There are direct dependencies of the global category on external packages, such as CLHEP, STL, and miscellaneous system utilities.

Within the management sub-category are “utility” classes generally used within the Geant4 kernel. They are, for the most part, uncorrelated with one another and include:

  • G4Allocator

  • G4FastVector

  • G4ReferenceCountedHandle

  • G4PhysicsVector, G4LPhysicsFreeVector, G4PhysicsOrderedFreeVector

  • G4Timer

  • G4UserLimits

  • G4UnitsTable

A general description of these classes is given in section 3.2 of the Geant4 User’s Guide for Application Developers.

The module includes wrappers to most CLHEP classes used in Geant4, and tools for memory management (G4Cache, G4AutoDelete) and for threading (G4AutoLock, G4Threading, G4ThreadLocalSingleton, G4TWorkspacePool). It also provides specialised fast implementations for some heavily used mathematical functions, like G4Exp, G4Log, G4Pow.

In applications where it is necessary to generate random numbers (normally from the same engine) in many different methods and parts of the program, it is highly desirable not to rely on or require knowledge of the global objects instantiated. By using static methods via a unique generator, the randomness of a sequence of numbers is best assured. Hence the use of a static generator has been introduced in the original design of HEPRandom as a project requirement in Geant4.

Class Design

Analysis and design of the HEPRandom module have been achieved following the Booch Object-Oriented methodology. Some of the original design diagrams in Booch notation are reported below. The figure is a general picture of the static class diagram.

  • HepRandomEngine - abstract class defining the interface for each Random engine. Its pure virtual methods must be defined by its subclasses representing the concrete Random engines.

  • HepJamesRandom - class inheriting from HepRandomEngine and defining a flat random number generator according to the algorithm described in “F.James, Comp.Phys.Comm. 60 (1990) 329”.

  • DRand48Engine - class inheriting from HepRandomEngine and defining a flat random number generator according to the drand48() and srand48() system functions from the C standard library.

  • RandEngine - class inheriting from HepRandomEngine and defining a flat random number generator according to the rand() and srand() system functions from the C standard library.

  • RanluxEngine - class inheriting from HepRandomEngine and defining a flat random number generator according to the algorithm described in “F.James, Comp.Phys.Comm. 60 (1990) 329-344” and originally implemented in FORTRAN 77 as part of the MATHLIB HEP library. It provides 5 different “luxury” levels [0..4].

  • RanecuEngine - class inheriting from HepRandomEngine and defining a flat random number generator according to the algorithm RANECU originally written in FORTRAN 77 as part of the MATHLIB HEP library. It uses a table of seeds which provides uncorrelated couples of seed values.

  • MixMaxRng - class inheriting from HepRandomEngine and interfacing the MixMax Matrix PseudoRandom Number Generator described in J.Comput.Phys. 97, 573 (1991). This class is instantiated by default as the default random engine.

  • HepRandom - the main class collecting all the methods defining the different random generators applied to HepRandomEngine. It is a singleton class which all the distribution classes derive from. This singleton is instantiated by default.

  • RandFlat - distribution class for flat random number generation. It also provides methods to fill an array of flat random values, given its size or shoot bits.

  • RandExponential - distribution class defining exponential random number distribution, given a mean. It also provides a method to fill an array of flat random values, given its size.

  • RandGauss - distribution class defining Gauss random number distribution, given a mean or specifying also a deviation. It also provides a method to fill an array of flat random values, given its size.

  • RandBreitWigner - distribution class defining the Breit-Wigner random number distribution. It also provides a method to fill an array of flat random values, given its size.

  • RandPoisson - distribution class defining Poisson random number distribution, given a mean. It also provides a method to fill an array of flat random values, given its size.

../../_images/classDgmRandom.jpg

Fig. 24 HEPRandom module

For detailed documentation about the HEPRandom classes see the CLHEP documentation <http://proj-clhep.web.cern.ch/proj-clhep/#docu>.

Information written in this manual is extracted from the original manifesto distributed with the HEPRandom package <http://proj-clhep.web.cern.ch/proj-clhep/doc/CLHEP_1_7/UserGuide/Random/Random.html>`__.

HEPNumerics

The HEPNumerics module includes a set of classes which implement numerical algorithms for general use in Geant4. The User’s Guide for Application Developers contains a description of each class. Most of the algorithms were implemented using methods from the following books:

  • B.H. Flowers, “An introduction to Numerical Methods In C++”, Clarendon Press, Oxford 1995.

  • M. Abramowitz, I. Stegun, “Handbook of mathematical functions”, DOVER Publications INC, New York 1965; chapters 9, 10, and 22.

The HEPNumerics module provides general mathematical methods supporting Geant4 Monte-Carlo simulation processes. Among these, there are methods for function and array interpolations using known special functions, class method integration solving polynomial equation (up to 4th order) and some others.

Of particular interest is the templated class G4Integrator which consists of methods allowing to integrate class methods. Since the type whose method should be integrated is not known in advance, G4Integrator uses templated signatures and pointers to functions in its API. It provides both usual numerical methods like adaptive Gauss or Simpson integration, and more sophisticated (faster and at the same accurate) methods based on the orthogonal polynomials.

Among the different integration methods involving orthogonal polynomials there are: Gauss-Legendre, Gauss-Chebyshev, Gauss-Hermite and Gauss-Jacobi methods:

template <class T, class F>
G4double G4Integrator<T,F>::Legendre( T& typeT, F f, G4double a,
                                    G4double b, G4int nLegendre )
//
// The value nLegendre set the accuracy required, i.e the number of points
// where the function pFunction will be evaluated during integration.
// The function creates the arrays for abscissas and weights that used
// in Gauss-Legendre quadrature method.
// The values a and b are the limits of integration of the function  f.
// nLegendre MUST BE EVEN !!!
// Returns the integral of the function f between a and b, by 2*fNumber point
// Gauss-Legendre integration: the function is evaluated exactly
// 2*fNumber times at interior points in the range of integration.
// Since the weights and abscissas are, in this case, symmetric around
// the midpoint of the range of integration, there are actually only
// fNumber distinct values of each.
// Convenient for using with some class object dataT

template <class T, class F>
G4double G4Integrator<T,F>::Legendre10( T& typeT, F f,G4double a,
                                    G4double b )
//
// Returns the integral of the function to be pointed by T::f between a and b,
// by ten point Gauss-Legendre integration: the function is evaluated exactly
// ten times at interior points in the range of integration. Since the weights
// and abscissas are, in this case, symmetric around the midpoint of the
// range of integration, there are actually only five distinct values of each
// Convenient for using with class object typeT. The method is very fast and accurate enough.
// Roots and weights are from Abramowitz M., Stegun I.A. 1964 , Handbook of Math... , p. 916

template <class T, class F>
G4double G4Integrator<T,F>::Legendre96( T& typeT, F f,G4double a,
                                    G4double b )
//
// Returns the integral of the function to be pointed by T::f between a and b,
// by 96 point Gauss-Legendre integration: the function is evaluated exactly
// ten Times at interior points in the range of integration. Since the weights
// and abscissas are, in this case, symmetric around the midpoint of the
// range of integration, there are actually only five distinct values of each
// Convenient for using with some class object typeT. The method is very accurate and fast enough.
// Roots and weights are from Abramowitz M., Stegun I.A. 1964 , Handbook of Math... , p. 919

template <class T, class F>
G4double G4Integrator<T,F>::Chebyshev( T& typeT, F f, G4double a,
                                    G4double b, G4int nChebyshev )
//
// Integrates function pointed by T::f from a to b by Gauss-Chebyshev
// quadrature method.
// Convenient for using with class object typeT

template <class T, class F>
G4double G4Integrator<T,F>::Laguerre( T& typeT, F f, G4double alpha,
                                    G4int nLaguerre )
//
// Integral from zero to infinity of std::pow(x,alpha)*std::exp(-x)*f(x).
// The value of nLaguerre sets the accuracy.
// The function creates arrays fAbscissa[0,..,nLaguerre-1] and
// fWeight[0,..,nLaguerre-1] .
// Convenient for using with class object 'typeT' and (typeT.*f) function
// (T::f)

template <class T, class F>
G4double G4Integrator<T,F>::Hermite( T& typeT, F f, G4int nHermite )
//
// Gauss-Hermite method for integration of std::exp(-x*x)*f(x)
// from minus infinity to plus infinity.

template <class T, class F>
G4double G4Integrator<T,F>::Jacobi( T& typeT, F f, G4double alpha,
                                    G4double beta, G4int nJacobi )
//
// Gauss-Jacobi method for integration of ((1-x)^alpha)*((1+x)^beta)*f(x)
// from minus unit to plus unit (-1,+1).

HEPGeometry

Documentation for the HEPGeometry module is provided in the CLHEP documentation.