Generate momenta using the RauboldLynch method. More...

#include <G4INCLPhaseSpaceRauboldLynch.hh>

Inheritance diagram for G4INCL::PhaseSpaceRauboldLynch:

Public Member Functions
void	generate (const G4double sqrtS, ParticleList &particles)
	Generate momenta according to a uniform, Lorentz-invariant phase-space model. More...

G4double	getMaxGeneratedWeight () const
	Return the largest generated weight. More...

PhaseSpaceRauboldLynch &	operator= (PhaseSpaceRauboldLynch const &rhs)
	Dummy assignment operator to silence Coverity warning. More...

	PhaseSpaceRauboldLynch ()

	PhaseSpaceRauboldLynch (PhaseSpaceRauboldLynch const &other)
	Dummy copy constructor to silence Coverity warning. More...

virtual	~PhaseSpaceRauboldLynch ()

Private Member Functions
G4double	computeMaximumWeightNaive ()
	Compute the maximum possible weight using a naive algorithm. More...

G4double	computeMaximumWeightParam ()
	Compute the maximum possible weight using parametrizations. More...

G4double	computeWeight ()
	Compute the maximum possible weight. More...

void	generateEvent (ParticleList &particles)
	Generate an event. More...

void	initialize (ParticleList &particles)
	Initialize internal structures (masses and sum of masses) More...

Private Attributes
G4double	availableEnergy

std::vector< G4double >	invariantMasses

std::vector< G4double >	masses

G4double	maxGeneratedWeight

std::vector< G4double >	momentaCM

size_t	nParticles

G4double	prelog [wMaxNP]
	Precalculated coefficients: -ln(n) More...

std::vector< G4double >	rnd

G4double	sqrtS

std::vector< G4double >	sumMasses

InterpolationTable *	wMaxCorrection

InterpolationTable *	wMaxMassless

Static Private Attributes
static const size_t	nMasslessParticlesTable = 13

static const G4double	wMaxCorrectionX [wMaxNE]

static const G4double	wMaxCorrectionY [wMaxNE]

static const G4double	wMaxInterpolationMargin = std::log(1.5)

static const G4double	wMaxMasslessX [wMaxNE]

static const G4double	wMaxMasslessY [wMaxNE]

static const size_t	wMaxNE = 30

static const size_t	wMaxNP = 20

Detailed Description

Generate momenta using the RauboldLynch method.

Definition at line 49 of file G4INCLPhaseSpaceRauboldLynch.hh.

Constructor & Destructor Documentation

◆ PhaseSpaceRauboldLynch() [1/2]

G4INCL::PhaseSpaceRauboldLynch::PhaseSpaceRauboldLynch ( )

Definition at line 180 of file G4INCLPhaseSpaceRauboldLynch.cc.

                                                 :
    nParticles(0),
    sqrtS(0.),
    availableEnergy(0.),
    maxGeneratedWeight(0.)
  {
    std::vector<G4double> wMaxMasslessXV(wMaxMasslessX, wMaxMasslessX + wMaxNE);
    std::vector<G4double> wMaxMasslessYV(wMaxMasslessY, wMaxMasslessY + wMaxNE);
    wMaxMassless = new InterpolationTable(wMaxMasslessXV, wMaxMasslessYV);
    std::vector<G4double> wMaxCorrectionXV(wMaxCorrectionX, wMaxCorrectionX + wMaxNE);
    std::vector<G4double> wMaxCorrectionYV(wMaxCorrectionY, wMaxCorrectionY + wMaxNE);
    wMaxCorrection = new InterpolationTable(wMaxCorrectionXV, wMaxCorrectionYV);
 
    // initialize the precalculated coefficients
    prelog[0] = 0.;
    for(size_t i=1; i<wMaxNP; ++i) {
      prelog[i] = -std::log(G4double(i));
    }
  }

References prelog, wMaxCorrection, wMaxCorrectionX, wMaxCorrectionY, wMaxMassless, wMaxMasslessX, wMaxMasslessY, wMaxNE, and wMaxNP.

◆ ~PhaseSpaceRauboldLynch()

G4INCL::PhaseSpaceRauboldLynch::~PhaseSpaceRauboldLynch ( )

virtual

Definition at line 200 of file G4INCLPhaseSpaceRauboldLynch.cc.

                                                  {
    delete wMaxMassless;
    delete wMaxCorrection;
  }

References wMaxCorrection, and wMaxMassless.

◆ PhaseSpaceRauboldLynch() [2/2]

G4INCL::PhaseSpaceRauboldLynch::PhaseSpaceRauboldLynch ( PhaseSpaceRauboldLynch const & other )

Dummy copy constructor to silence Coverity warning.

Member Function Documentation

◆ computeMaximumWeightNaive()

G4double G4INCL::PhaseSpaceRauboldLynch::computeMaximumWeightNaive ( )

private

Compute the maximum possible weight using a naive algorithm.

Definition at line 257 of file G4INCLPhaseSpaceRauboldLynch.cc.

                                                             {
    // compute the maximum weight
    G4double eMMax = sqrtS + masses[0];
    G4double eMMin = 0.;
    G4double wMax = 1.;
    for(size_t i=1; i<nParticles; i++) {
      eMMin += masses[i-1];
      eMMax += masses[i];
      wMax *= KinematicsUtils::momentumInCM(eMMax, eMMin, masses[i]);
    }
    return wMax;
  }

References masses, G4INCL::KinematicsUtils::momentumInCM(), nParticles, and sqrtS.

Referenced by computeMaximumWeightParam().

◆ computeMaximumWeightParam()

G4double G4INCL::PhaseSpaceRauboldLynch::computeMaximumWeightParam ( )

private

Compute the maximum possible weight using parametrizations.

Definition at line 270 of file G4INCLPhaseSpaceRauboldLynch.cc.

                                                             {
#ifndef INCLXX_IN_GEANT4_MODE
    if(nParticles>=wMaxNP) {
      INCL_WARN("The requested number of particles (" << nParticles << ") requires extrapolation the tables in PhaseSpaceRauboldLynch. YMMV." << '\n');
    }
 
    if(availableEnergy>wMaxMasslessX[wMaxNE-1] || availableEnergy<wMaxMasslessX[0] ) {
      INCL_WARN("The requested available energy (" << availableEnergy << " MeV) requires extrapolation the tables in PhaseSpaceRauboldLynch. YMMV." << '\n');
    }
#endif
 
    // compute the maximum weight
    const G4double logMassless = ((*wMaxMassless)(availableEnergy)+prelog[nParticles])*(nParticles-1);
    const G4double reducedSqrtS = availableEnergy/sumMasses[nParticles-1];
    const G4double correction = (*wMaxCorrection)(reducedSqrtS);
    const G4double wMax = std::exp(correction*G4double(nParticles-1) + logMassless + wMaxInterpolationMargin);
    if(wMax>0.)
      return wMax;
    else
      return computeMaximumWeightNaive();
  }

References availableEnergy, computeMaximumWeightNaive(), INCL_WARN, nParticles, prelog, sumMasses, wMaxInterpolationMargin, wMaxMasslessX, wMaxNE, and wMaxNP.

Referenced by generate().

◆ computeWeight()

G4double G4INCL::PhaseSpaceRauboldLynch::computeWeight ( )

private

Compute the maximum possible weight.

Definition at line 292 of file G4INCLPhaseSpaceRauboldLynch.cc.

                                                 {
    // Generate nParticles-2 sorted random numbers, bracketed by 0 and 1
    rnd[0] = 0.;
    for(size_t i=1; i<nParticles-1; ++i)
      rnd[i] = Random::shoot();
    rnd[nParticles-1] = 1.;
    std::sort(rnd.begin()+1, rnd.begin()+nParticles-1);
    
    // invariant masses
    for(size_t i=0; i<nParticles; ++i)
      invariantMasses[i] = rnd[i]*availableEnergy + sumMasses[i];
 
    // compute the CM momenta and the weight for the current event
    G4double weight=KinematicsUtils::momentumInCM(invariantMasses[1], invariantMasses[0], masses[1]);
    momentaCM[0] = weight;
    for(size_t i=1; i<nParticles-1; ++i) {
      G4double momentumCM;
// assert(invariantMasses[i+1]-invariantMasses[i]-masses[i+1]> - 1.E-8);
      if(invariantMasses[i+1]-invariantMasses[i]-masses[i+1] < 0.) momentumCM = 0.;
      else momentumCM = KinematicsUtils::momentumInCM(invariantMasses[i+1], invariantMasses[i], masses[i+1]);
      momentaCM[i] = momentumCM;
      weight *= momentumCM;
    }
 
    return weight;
  }

References availableEnergy, invariantMasses, masses, momentaCM, G4INCL::KinematicsUtils::momentumInCM(), nParticles, rnd, G4INCL::Random::shoot(), and sumMasses.

Referenced by generate().

◆ generate()

void G4INCL::PhaseSpaceRauboldLynch::generate	(	const G4double	sqrtS,
		ParticleList &	particles
	)

virtual

Generate momenta according to a uniform, Lorentz-invariant phase-space model.

This function will assign momenta to the particles in the list that is passed as an argument. The event is generated in the CM frame.

Parameters

sqrtS	total centre-of-mass energy of the system
particles	list of particles

Implements G4INCL::IPhaseSpaceGenerator.

Definition at line 205 of file G4INCLPhaseSpaceRauboldLynch.cc.

                                                                                   {
    maxGeneratedWeight = 0.;
 
    sqrtS = sqs;
 
    // initialize the structures containing the particle masses
    initialize(particles);
 
    // initialize the maximum weight
    const G4double weightMax = computeMaximumWeightParam();
 
    const G4int maxIter = 500;
    G4int iter = 0;
    G4double weight, r;
    do {
      weight = computeWeight();
      maxGeneratedWeight = std::max(weight, maxGeneratedWeight);
// assert(weight<=weightMax);
      r = Random::shoot();
    } while(++iter<maxIter && r*weightMax>weight); /* Loop checking, 10.07.2015, D.Mancusi */
 
#ifndef INCLXX_IN_GEANT4_MODE
    if(iter>=maxIter) {
      INCL_WARN("Number of tries exceeded in PhaseSpaceRauboldLynch::generate()\n"
                << "nParticles=" << nParticles << ", sqrtS=" << sqrtS << '\n');
    }
#endif
 
    generateEvent(particles);
  }

References computeMaximumWeightParam(), computeWeight(), generateEvent(), INCL_WARN, initialize(), G4INCL::Math::max(), maxGeneratedWeight, nParticles, G4INCL::Random::shoot(), and sqrtS.

◆ generateEvent()

void G4INCL::PhaseSpaceRauboldLynch::generateEvent ( ParticleList & particles )

private

Generate an event.

Definition at line 319 of file G4INCLPhaseSpaceRauboldLynch.cc.

                                                                    {
    Particle *p = particles[0];
    ThreeVector momentum = Random::normVector(momentaCM[0]);
    p->setMomentum(momentum);
    p->adjustEnergyFromMomentum();
 
    ThreeVector boostV;
 
    for(size_t i=1; i<nParticles; ++i) {
      p = particles[i];
      p->setMomentum(-momentum);
      p->adjustEnergyFromMomentum();
 
      if(i==nParticles-1)
        break;
 
      momentum = Random::normVector(momentaCM[i]);
 
      const G4double iM = invariantMasses[i];
      const G4double recoilE = std::sqrt(momentum.mag2() + iM*iM);
      boostV = -momentum/recoilE;
      for(size_t j=0; j<=i; ++j)
        particles[j]->boost(boostV);
    }
  }

References G4INCL::Particle::adjustEnergyFromMomentum(), invariantMasses, G4INCL::ThreeVector::mag2(), momentaCM, G4INCL::Random::normVector(), nParticles, and G4INCL::Particle::setMomentum().

Referenced by generate().

◆ getMaxGeneratedWeight()

G4double G4INCL::PhaseSpaceRauboldLynch::getMaxGeneratedWeight ( ) const

Return the largest generated weight.

Definition at line 345 of file G4INCLPhaseSpaceRauboldLynch.cc.

                                                               {
    return maxGeneratedWeight;
  }

References maxGeneratedWeight.

◆ initialize()

void G4INCL::PhaseSpaceRauboldLynch::initialize ( ParticleList & particles )

private

Initialize internal structures (masses and sum of masses)

Definition at line 236 of file G4INCLPhaseSpaceRauboldLynch.cc.

                                                                 {
    nParticles = particles.size();
// assert(nParticles>2);
 
    // masses and sum of masses
    masses.resize(nParticles);
    sumMasses.resize(nParticles);
    std::transform(particles.begin(), particles.end(), masses.begin(), std::mem_fn(&Particle::getMass));
    std::partial_sum(masses.begin(), masses.end(), sumMasses.begin());
 
    // sanity check
    availableEnergy  = sqrtS-sumMasses[nParticles-1];
// assert(availableEnergy>-1.e-5);
    if(availableEnergy<0.)
      availableEnergy = 0.;
 
    rnd.resize(nParticles);
    invariantMasses.resize(nParticles);
    momentaCM.resize(nParticles-1);
  }

References availableEnergy, G4INCL::Particle::getMass(), invariantMasses, masses, momentaCM, nParticles, rnd, sqrtS, sumMasses, and G4coutFormatters::anonymous_namespace{G4coutFormatters.cc}::transform().

Referenced by generate().

◆ operator=()

PhaseSpaceRauboldLynch & G4INCL::PhaseSpaceRauboldLynch::operator= ( PhaseSpaceRauboldLynch const & rhs )

Dummy assignment operator to silence Coverity warning.

Field Documentation

◆ availableEnergy

G4double G4INCL::PhaseSpaceRauboldLynch::availableEnergy

private

Definition at line 82 of file G4INCLPhaseSpaceRauboldLynch.hh.

Referenced by computeMaximumWeightParam(), computeWeight(), and initialize().

◆ invariantMasses

std::vector<G4double> G4INCL::PhaseSpaceRauboldLynch::invariantMasses

private

Definition at line 78 of file G4INCLPhaseSpaceRauboldLynch.hh.

Referenced by computeWeight(), generateEvent(), and initialize().

◆ masses

std::vector<G4double> G4INCL::PhaseSpaceRauboldLynch::masses

private

Definition at line 75 of file G4INCLPhaseSpaceRauboldLynch.hh.

Referenced by computeMaximumWeightNaive(), computeWeight(), and initialize().

◆ maxGeneratedWeight

G4double G4INCL::PhaseSpaceRauboldLynch::maxGeneratedWeight

private

Definition at line 83 of file G4INCLPhaseSpaceRauboldLynch.hh.

Referenced by generate(), and getMaxGeneratedWeight().

◆ momentaCM

std::vector<G4double> G4INCL::PhaseSpaceRauboldLynch::momentaCM

private

Definition at line 79 of file G4INCLPhaseSpaceRauboldLynch.hh.

Referenced by computeWeight(), generateEvent(), and initialize().

◆ nMasslessParticlesTable

const size_t G4INCL::PhaseSpaceRauboldLynch::nMasslessParticlesTable = 13

staticprivate

Definition at line 85 of file G4INCLPhaseSpaceRauboldLynch.hh.

◆ nParticles

size_t G4INCL::PhaseSpaceRauboldLynch::nParticles

private

Definition at line 80 of file G4INCLPhaseSpaceRauboldLynch.hh.

Referenced by computeMaximumWeightNaive(), computeMaximumWeightParam(), computeWeight(), generate(), generateEvent(), and initialize().

◆ prelog

G4double G4INCL::PhaseSpaceRauboldLynch::prelog[wMaxNP]

private

Precalculated coefficients: -ln(n)

Definition at line 99 of file G4INCLPhaseSpaceRauboldLynch.hh.

Referenced by computeMaximumWeightParam(), and PhaseSpaceRauboldLynch().

◆ rnd

std::vector<G4double> G4INCL::PhaseSpaceRauboldLynch::rnd

private

Definition at line 77 of file G4INCLPhaseSpaceRauboldLynch.hh.

Referenced by computeWeight(), and initialize().

◆ sqrtS

G4double G4INCL::PhaseSpaceRauboldLynch::sqrtS

private

Definition at line 81 of file G4INCLPhaseSpaceRauboldLynch.hh.

Referenced by computeMaximumWeightNaive(), generate(), and initialize().

◆ sumMasses

std::vector<G4double> G4INCL::PhaseSpaceRauboldLynch::sumMasses

private

Definition at line 76 of file G4INCLPhaseSpaceRauboldLynch.hh.

Referenced by computeMaximumWeightParam(), computeWeight(), and initialize().

◆ wMaxCorrection

InterpolationTable* G4INCL::PhaseSpaceRauboldLynch::wMaxCorrection

private

Definition at line 94 of file G4INCLPhaseSpaceRauboldLynch.hh.

Referenced by PhaseSpaceRauboldLynch(), and ~PhaseSpaceRauboldLynch().

◆ wMaxCorrectionX

const G4double G4INCL::PhaseSpaceRauboldLynch::wMaxCorrectionX

staticprivate

Definition at line 89 of file G4INCLPhaseSpaceRauboldLynch.hh.

Referenced by PhaseSpaceRauboldLynch().

◆ wMaxCorrectionY

const G4double G4INCL::PhaseSpaceRauboldLynch::wMaxCorrectionY

staticprivate

Definition at line 90 of file G4INCLPhaseSpaceRauboldLynch.hh.

Referenced by PhaseSpaceRauboldLynch().

◆ wMaxInterpolationMargin

const G4double G4INCL::PhaseSpaceRauboldLynch::wMaxInterpolationMargin = std::log(1.5)

staticprivate

Definition at line 91 of file G4INCLPhaseSpaceRauboldLynch.hh.

Referenced by computeMaximumWeightParam().

◆ wMaxMassless

InterpolationTable* G4INCL::PhaseSpaceRauboldLynch::wMaxMassless

private

Definition at line 93 of file G4INCLPhaseSpaceRauboldLynch.hh.

Referenced by PhaseSpaceRauboldLynch(), and ~PhaseSpaceRauboldLynch().

◆ wMaxMasslessX

const G4double G4INCL::PhaseSpaceRauboldLynch::wMaxMasslessX

staticprivate

Definition at line 87 of file G4INCLPhaseSpaceRauboldLynch.hh.

Referenced by computeMaximumWeightParam(), and PhaseSpaceRauboldLynch().

◆ wMaxMasslessY

const G4double G4INCL::PhaseSpaceRauboldLynch::wMaxMasslessY

staticprivate

Definition at line 88 of file G4INCLPhaseSpaceRauboldLynch.hh.

Referenced by PhaseSpaceRauboldLynch().

◆ wMaxNE

const size_t G4INCL::PhaseSpaceRauboldLynch::wMaxNE = 30

staticprivate

Definition at line 86 of file G4INCLPhaseSpaceRauboldLynch.hh.

Referenced by computeMaximumWeightParam(), and PhaseSpaceRauboldLynch().

◆ wMaxNP

const size_t G4INCL::PhaseSpaceRauboldLynch::wMaxNP = 20

staticprivate

Definition at line 96 of file G4INCLPhaseSpaceRauboldLynch.hh.

Referenced by computeMaximumWeightParam(), and PhaseSpaceRauboldLynch().

The documentation for this class was generated from the following files:

source/processes/hadronic/models/inclxx/incl_physics/include/G4INCLPhaseSpaceRauboldLynch.hh
source/processes/hadronic/models/inclxx/incl_physics/src/G4INCLPhaseSpaceRauboldLynch.cc

Public Member Functions

Private Member Functions

Private Attributes

Static Private Attributes

Detailed Description

Constructor & Destructor Documentation

◆ PhaseSpaceRauboldLynch() [1/2]

◆ ~PhaseSpaceRauboldLynch()

◆ PhaseSpaceRauboldLynch() [2/2]

Member Function Documentation

◆ computeMaximumWeightNaive()

◆ computeMaximumWeightParam()

◆ computeWeight()

◆ generate()

◆ generateEvent()

◆ getMaxGeneratedWeight()

◆ initialize()

◆ operator=()

Field Documentation

◆ availableEnergy

◆ invariantMasses

◆ masses

◆ maxGeneratedWeight

◆ momentaCM

◆ nMasslessParticlesTable

◆ nParticles

◆ prelog

◆ rnd

◆ sqrtS

◆ sumMasses

◆ wMaxCorrection

◆ wMaxCorrectionX

◆ wMaxCorrectionY

◆ wMaxInterpolationMargin

◆ wMaxMassless

◆ wMaxMasslessX

◆ wMaxMasslessY

◆ wMaxNE

◆ wMaxNP