Geant4-11
Public Member Functions | Protected Member Functions | Protected Attributes | Private Member Functions | Private Attributes
G4HETCAlpha Class Reference

#include <G4HETCAlpha.hh>

Inheritance diagram for G4HETCAlpha:
G4HETCChargedFragment G4HETCFragment G4VPreCompoundFragment

Public Member Functions

G4double CalcEmissionProbability (const G4Fragment &aFragment)
 
 G4HETCAlpha ()
 
G4int GetA () const
 
G4double GetBindingEnergy () const
 
G4double GetEmissionProbability () const
 
G4double GetEnergyThreshold () const
 
const G4LorentzVectorGetMomentum () const
 
G4double GetNuclearMass () const
 
G4ReactionProductGetReactionProduct () const
 
G4int GetRestA () const
 
G4double GetRestNuclearMass () const
 
G4int GetRestZ () const
 
G4int GetZ () const
 
void Initialize (const G4Fragment &aFragment)
 
G4bool IsItPossible (const G4Fragment &aFragment) const
 
G4bool operator!= (const G4VPreCompoundFragment &right) const =delete
 
G4bool operator== (const G4VPreCompoundFragment &right) const =delete
 
virtual G4double SampleKineticEnergy (const G4Fragment &aFragment)
 
void SetMomentum (const G4LorentzVector &value)
 
void SetOPTxs (G4int)
 
void UseSICB (G4bool)
 
 ~G4HETCAlpha ()
 

Protected Member Functions

G4double BetaRand (G4int N, G4int L) const
 
virtual G4double GetAlpha () const
 
virtual G4double GetBeta () const
 
virtual G4double GetSpinFactor () const
 
virtual G4double K (const G4Fragment &aFragment)
 

Protected Attributes

G4NuclearLevelDatafNucData
 
G4Powg4calc
 
G4int OPTxs
 
G4int theA
 
G4double theBindingEnergy
 
G4double theCoulombBarrier
 
G4double theEmissionProbability
 
G4int theFragA
 
G4int theFragZ
 
G4double theMass
 
G4double theMaxKinEnergy
 
G4double theMinKinEnergy
 
G4DeexPrecoParameterstheParameters
 
G4double theReducedMass
 
G4int theResA
 
G4double theResA13
 
G4double theResMass
 
G4int theResZ
 
G4int theZ
 
G4bool useSICB
 

Private Member Functions

 G4HETCAlpha (const G4HETCAlpha &right)
 
G4double IntegrateEmissionProbability (G4double &Low, G4double &Up, const G4Fragment &aFragment)
 
G4bool operator!= (const G4HETCAlpha &right) const
 
const G4HETCAlphaoperator= (const G4HETCAlpha &right)
 
G4bool operator== (const G4HETCAlpha &right) const
 

Private Attributes

const G4ParticleDefinitionparticle
 
G4double r2norm
 
G4AlphaCoulombBarrier theAlphaCoulombBarrier
 
G4VCoulombBarriertheCoulombBarrierPtr
 
G4LorentzVector theMomentum
 

Detailed Description

Definition at line 40 of file G4HETCAlpha.hh.

Constructor & Destructor Documentation

◆ G4HETCAlpha() [1/2]

G4HETCAlpha::G4HETCAlpha ( )

Definition at line 36 of file G4HETCAlpha.cc.

38{}
static G4Alpha * Alpha()
Definition: G4Alpha.cc:88
G4AlphaCoulombBarrier theAlphaCoulombBarrier
Definition: G4HETCAlpha.hh:66

◆ ~G4HETCAlpha()

G4HETCAlpha::~G4HETCAlpha ( )

Definition at line 40 of file G4HETCAlpha.cc.

41{}

◆ G4HETCAlpha() [2/2]

G4HETCAlpha::G4HETCAlpha ( const G4HETCAlpha right)
private

Member Function Documentation

◆ BetaRand()

G4double G4HETCFragment::BetaRand ( G4int  N,
G4int  L 
) const
inlineprotectedinherited

Definition at line 78 of file G4HETCFragment.hh.

79{
82
83 return Y1/(Y1+Y2);
84}
static constexpr double L
Definition: G4SIunits.hh:104
double G4double
Definition: G4Types.hh:83
ThreeVector shoot(const G4int Ap, const G4int Af)

References L, and G4INCL::DeJongSpin::shoot().

Referenced by G4HETCChargedFragment::SampleKineticEnergy(), and G4HETCNeutron::SampleKineticEnergy().

◆ CalcEmissionProbability()

G4double G4HETCFragment::CalcEmissionProbability ( const G4Fragment aFragment)
virtualinherited

Implements G4VPreCompoundFragment.

Definition at line 50 of file G4HETCFragment.cc.

52{
53 if (GetEnergyThreshold() <= 0.0)
54 {
56 return 0.0;
57 }
58 // Coulomb barrier is the lower limit
59 // of integration over kinetic energy
62
64}
G4double IntegrateEmissionProbability(G4double &Low, G4double &Up, const G4Fragment &aFragment)
G4double GetEnergyThreshold() const

References G4VPreCompoundFragment::GetEnergyThreshold(), G4HETCFragment::IntegrateEmissionProbability(), G4VPreCompoundFragment::theCoulombBarrier, G4VPreCompoundFragment::theEmissionProbability, and G4VPreCompoundFragment::theMaxKinEnergy.

◆ GetA()

G4int G4VPreCompoundFragment::GetA ( ) const
inlineinherited

◆ GetAlpha()

G4double G4HETCAlpha::GetAlpha ( ) const
protectedvirtual

Implements G4HETCFragment.

Definition at line 43 of file G4HETCAlpha.cc.

44{
45 G4double C = 0.0;
46 if (theFragZ <= 30)
47 {
48 C = 0.10;
49 }
50 else if (theFragZ <= 50)
51 {
52 C = 0.1 + -((theFragZ-50.)/20.)*0.02;
53 }
54 else if (theFragZ < 70)
55 {
56 C = 0.08 + -((theFragZ-70.)/20.)*0.02;
57 }
58 else
59 {
60 C = 0.06;
61 }
62 return 1.0+C;
63}
G4double C(G4double temp)

References C(), and G4VPreCompoundFragment::theFragZ.

◆ GetBeta()

G4double G4HETCAlpha::GetBeta ( ) const
protectedvirtual

Implements G4HETCFragment.

Definition at line 65 of file G4HETCAlpha.cc.

66{
67 return theCoulombBarrier;
68}

References G4VPreCompoundFragment::theCoulombBarrier.

◆ GetBindingEnergy()

G4double G4VPreCompoundFragment::GetBindingEnergy ( ) const
inlineinherited

◆ GetEmissionProbability()

G4double G4VPreCompoundFragment::GetEmissionProbability ( ) const
inlineinherited

◆ GetEnergyThreshold()

G4double G4VPreCompoundFragment::GetEnergyThreshold ( ) const
inlineinherited

◆ GetMomentum()

const G4LorentzVector & G4VPreCompoundFragment::GetMomentum ( ) const
inlineinherited

◆ GetNuclearMass()

G4double G4VPreCompoundFragment::GetNuclearMass ( ) const
inlineinherited

◆ GetReactionProduct()

G4ReactionProduct * G4VPreCompoundFragment::GetReactionProduct ( ) const
inlineinherited

◆ GetRestA()

G4int G4VPreCompoundFragment::GetRestA ( ) const
inlineinherited

◆ GetRestNuclearMass()

G4double G4VPreCompoundFragment::GetRestNuclearMass ( ) const
inlineinherited

◆ GetRestZ()

G4int G4VPreCompoundFragment::GetRestZ ( ) const
inlineinherited

◆ GetSpinFactor()

G4double G4HETCAlpha::GetSpinFactor ( ) const
protectedvirtual

Implements G4HETCFragment.

Definition at line 70 of file G4HETCAlpha.cc.

71{
72 return 1.0;
73}

◆ GetZ()

G4int G4VPreCompoundFragment::GetZ ( ) const
inlineinherited

◆ Initialize()

void G4VPreCompoundFragment::Initialize ( const G4Fragment aFragment)
inherited

Definition at line 79 of file G4VPreCompoundFragment.cc.

80{
81 theFragA = aFragment.GetA_asInt();
82 theFragZ = aFragment.GetZ_asInt();
85
87 if ((theResA < theResZ) || (theResA < theA) || (theResZ < theZ)) {
88 return;
89 }
90
93 GetCoulombBarrier(theResA,theResZ,aFragment.GetExcitationEnergy());
94
96
97 // Calculate masses
100
101 // Compute Binding Energies for fragments
102 // needed to separate a fragment from the nucleus
104
105 // Compute Maximal Kinetic Energy which can be carried by fragments
106 // after separation - the true assimptotic value
107 G4double Ecm = aFragment.GetMomentum().m();
108 G4double twoEcm = Ecm + Ecm;
110 /twoEcm - theMass,0.0);
111 theMinKinEnergy = (elim == 0.0) ? 0.0 :
112 std::max(((theMass+elim)*(twoEcm-theMass-elim) +
113 theMass*theMass)/twoEcm - theMass,0.0);
114}
G4double GetGroundStateMass() const
Definition: G4Fragment.hh:304
G4double GetExcitationEnergy() const
Definition: G4Fragment.hh:299
const G4LorentzVector & GetMomentum() const
Definition: G4Fragment.hh:323
G4int GetZ_asInt() const
Definition: G4Fragment.hh:276
G4int GetA_asInt() const
Definition: G4Fragment.hh:271
static G4double GetNuclearMass(const G4double A, const G4double Z)
G4double Z13(G4int Z) const
Definition: G4Pow.hh:123
G4VCoulombBarrier * theCoulombBarrierPtr
T max(const T t1, const T t2)
brief Return the largest of the two arguments

References G4VPreCompoundFragment::g4calc, G4Fragment::GetA_asInt(), G4Fragment::GetExcitationEnergy(), G4Fragment::GetGroundStateMass(), G4Fragment::GetMomentum(), G4NucleiProperties::GetNuclearMass(), G4Fragment::GetZ_asInt(), CLHEP::HepLorentzVector::m(), G4INCL::Math::max(), G4VPreCompoundFragment::OPTxs, G4VPreCompoundFragment::theA, G4VPreCompoundFragment::theBindingEnergy, G4VPreCompoundFragment::theCoulombBarrier, G4VPreCompoundFragment::theCoulombBarrierPtr, G4VPreCompoundFragment::theFragA, G4VPreCompoundFragment::theFragZ, G4VPreCompoundFragment::theMass, G4VPreCompoundFragment::theMaxKinEnergy, G4VPreCompoundFragment::theMinKinEnergy, G4VPreCompoundFragment::theReducedMass, G4VPreCompoundFragment::theResA, G4VPreCompoundFragment::theResA13, G4VPreCompoundFragment::theResMass, G4VPreCompoundFragment::theResZ, G4VPreCompoundFragment::theZ, and G4Pow::Z13().

◆ IntegrateEmissionProbability()

G4double G4HETCFragment::IntegrateEmissionProbability ( G4double Low,
G4double Up,
const G4Fragment aFragment 
)
privateinherited

Definition at line 66 of file G4HETCFragment.cc.

69{
70 G4double U = aFragment.GetExcitationEnergy();
71
72 G4int P = aFragment.GetNumberOfParticles();
73 G4int H = aFragment.GetNumberOfHoles();
74 G4int N = P + H;
75 G4int Pb = P - theA;
76 G4int Nb = Pb + H;
77 if (Nb <= 0.0) { return 0.0; }
78
81
82 G4double A = G4double(P*P+H*H+P-3*H)/(4.0*ga);
83 G4double Ab = G4double(Pb*Pb+H*H+Pb-3*H)/(4.0*gb);
84 U = std::max(U-A,0.0);
85 if (U <= 0.0) { return 0.0; }
86
87 G4int Pf = P;
88 G4int Hf = H;
89 G4int Nf = N-1;
90 for (G4int i = 1; i < theA; ++i)
91 {
92 Pf *= (P-i);
93 Hf *= (H-i);
94 Nf *= (N-1-i);
95 }
96
97 G4double X = std::max(Up - Ab + GetBeta(),0.0);
98 G4double Y = std::max(Up - Ab - Low, 0.0);
99
101 *g4calc->Z23(theResA)*Pf*Hf*Nf*K(aFragment)*(X/Nb - Y/(Nb+1))
102 *U*g4calc->powN(gb*Y,Nb)/g4calc->powN(ga*U,N);
103
104 return Probability;
105}
G4double Y(G4double density)
static constexpr double pi2
Definition: G4SIunits.hh:58
int G4int
Definition: G4Types.hh:85
const G4double A[17]
G4int GetNumberOfParticles() const
Definition: G4Fragment.hh:361
G4int GetNumberOfHoles() const
Definition: G4Fragment.hh:381
virtual G4double GetSpinFactor() const =0
virtual G4double K(const G4Fragment &aFragment)=0
virtual G4double GetAlpha() const =0
virtual G4double GetBeta() const =0
G4double GetLevelDensity(G4int Z, G4int A, G4double U)
G4double powN(G4double x, G4int n) const
Definition: G4Pow.cc:166
G4double Z23(G4int Z) const
Definition: G4Pow.hh:125
G4NuclearLevelData * fNucData
const G4double Pf
Fermi momentum [MeV/c].
static double P[]

References A, G4VPreCompoundFragment::fNucData, G4VPreCompoundFragment::g4calc, G4HETCFragment::GetAlpha(), G4HETCFragment::GetBeta(), G4Fragment::GetExcitationEnergy(), G4NuclearLevelData::GetLevelDensity(), G4Fragment::GetNumberOfHoles(), G4Fragment::GetNumberOfParticles(), G4HETCFragment::GetSpinFactor(), G4HETCFragment::K(), G4INCL::Math::max(), P, G4INCL::PhysicalConstants::Pf, pi2, G4Pow::powN(), G4HETCFragment::r2norm, G4VPreCompoundFragment::theA, G4VPreCompoundFragment::theFragA, G4VPreCompoundFragment::theFragZ, G4VPreCompoundFragment::theReducedMass, G4VPreCompoundFragment::theResA, G4VPreCompoundFragment::theResZ, Y(), and G4Pow::Z23().

Referenced by G4HETCFragment::CalcEmissionProbability().

◆ IsItPossible()

G4bool G4VPreCompoundFragment::IsItPossible ( const G4Fragment aFragment) const
inlineinherited

◆ K()

G4double G4HETCAlpha::K ( const G4Fragment aFragment)
protectedvirtual

Implements G4HETCFragment.

Definition at line 75 of file G4HETCAlpha.cc.

76{
77 // Number of protons in emitted fragment
78 G4int Pa = theZ;
79 // Number of neutrons in emitted fragment
80 G4int Na = theA - Pa;
81
83
84 G4int P = aFragment.GetNumberOfParticles();
85 G4int H = aFragment.GetNumberOfHoles();
86
87 G4double result = 0.0;
88 if (P > 3)
89 {
90 result = 3.0/(P*(P-1.0)*(P-2.0)*(P-3.0)) *
91 (H*(H-1.0)*(H-2.0)*(H-3.0)*r*r*(r-1.0)*(r-1.0) +
92 2.0*H*(H-1.0)*(H-2.0)*(Pa*r*(1.0-r)*(1.0-r)+Na*r*r*(1.0-r)) +
93 H*(H-1.0)*(Pa*(Pa-1.0)*(1.0-r)*(1.0-r)+4.0*Na*Pa*r*(1.0-r)+Na*(Na-1.0)*r*r) +
94 2*H*(Pa*Na*(Na-1.0)*r+Pa*(Pa-1.0)*Na*(1.0-r)) +
95 Pa*(Pa-1.0)*Na*(Na-1.0));
96
97 result /= 6.0*r*r*(1. - r) *(1. - r);
98 }
99 return std::max(0.0,result);
100}

References G4Fragment::GetNumberOfHoles(), G4Fragment::GetNumberOfParticles(), G4INCL::Math::max(), P, G4VPreCompoundFragment::theA, G4VPreCompoundFragment::theResA, G4VPreCompoundFragment::theResZ, and G4VPreCompoundFragment::theZ.

◆ operator!=() [1/2]

G4bool G4HETCAlpha::operator!= ( const G4HETCAlpha right) const
private

◆ operator!=() [2/2]

G4bool G4VPreCompoundFragment::operator!= ( const G4VPreCompoundFragment right) const
deleteinherited

◆ operator=()

const G4HETCAlpha & G4HETCAlpha::operator= ( const G4HETCAlpha right)
private

◆ operator==() [1/2]

G4bool G4HETCAlpha::operator== ( const G4HETCAlpha right) const
private

◆ operator==() [2/2]

G4bool G4VPreCompoundFragment::operator== ( const G4VPreCompoundFragment right) const
deleteinherited

◆ SampleKineticEnergy()

G4double G4HETCChargedFragment::SampleKineticEnergy ( const G4Fragment aFragment)
virtualinherited

◆ SetMomentum()

void G4VPreCompoundFragment::SetMomentum ( const G4LorentzVector value)
inlineinherited

◆ SetOPTxs()

void G4VPreCompoundFragment::SetOPTxs ( G4int  )
inlineinherited

◆ UseSICB()

void G4VPreCompoundFragment::UseSICB ( G4bool  )
inlineinherited

Field Documentation

◆ fNucData

G4NuclearLevelData* G4VPreCompoundFragment::fNucData
protectedinherited

◆ g4calc

G4Pow* G4VPreCompoundFragment::g4calc
protectedinherited

◆ OPTxs

G4int G4VPreCompoundFragment::OPTxs
protectedinherited

◆ particle

const G4ParticleDefinition* G4VPreCompoundFragment::particle
privateinherited

◆ r2norm

G4double G4HETCFragment::r2norm
privateinherited

◆ theA

G4int G4VPreCompoundFragment::theA
protectedinherited

◆ theAlphaCoulombBarrier

G4AlphaCoulombBarrier G4HETCAlpha::theAlphaCoulombBarrier
private

Definition at line 66 of file G4HETCAlpha.hh.

◆ theBindingEnergy

G4double G4VPreCompoundFragment::theBindingEnergy
protectedinherited

◆ theCoulombBarrier

G4double G4VPreCompoundFragment::theCoulombBarrier
protectedinherited

◆ theCoulombBarrierPtr

G4VCoulombBarrier* G4VPreCompoundFragment::theCoulombBarrierPtr
privateinherited

Definition at line 132 of file G4VPreCompoundFragment.hh.

Referenced by G4VPreCompoundFragment::Initialize().

◆ theEmissionProbability

G4double G4VPreCompoundFragment::theEmissionProbability
protectedinherited

◆ theFragA

G4int G4VPreCompoundFragment::theFragA
protectedinherited

◆ theFragZ

G4int G4VPreCompoundFragment::theFragZ
protectedinherited

◆ theMass

G4double G4VPreCompoundFragment::theMass
protectedinherited

◆ theMaxKinEnergy

G4double G4VPreCompoundFragment::theMaxKinEnergy
protectedinherited

◆ theMinKinEnergy

G4double G4VPreCompoundFragment::theMinKinEnergy
protectedinherited

◆ theMomentum

G4LorentzVector G4VPreCompoundFragment::theMomentum
privateinherited

Definition at line 134 of file G4VPreCompoundFragment.hh.

◆ theParameters

G4DeexPrecoParameters* G4VPreCompoundFragment::theParameters
protectedinherited

◆ theReducedMass

G4double G4VPreCompoundFragment::theReducedMass
protectedinherited

◆ theResA

G4int G4VPreCompoundFragment::theResA
protectedinherited

◆ theResA13

G4double G4VPreCompoundFragment::theResA13
protectedinherited

◆ theResMass

G4double G4VPreCompoundFragment::theResMass
protectedinherited

Definition at line 153 of file G4VPreCompoundFragment.hh.

Referenced by G4VPreCompoundFragment::Initialize().

◆ theResZ

G4int G4VPreCompoundFragment::theResZ
protectedinherited

◆ theZ

G4int G4VPreCompoundFragment::theZ
protectedinherited

◆ useSICB

G4bool G4VPreCompoundFragment::useSICB
protectedinherited

Definition at line 163 of file G4VPreCompoundFragment.hh.


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