Geant4-11
Public Member Functions | Private Attributes
G4IonisParamElm Class Reference

#include <G4IonisParamElm.hh>

Public Member Functions

 G4IonisParamElm (__void__ &)
 
 G4IonisParamElm (G4double Z)
 
 G4IonisParamElm (G4IonisParamElm &)=delete
 
G4double GetAlow () const
 
G4double GetBlow () const
 
G4double GetClow () const
 
G4double GetFermiVelocity () const
 
G4double GetLFactor () const
 
G4double GetlogZ3 () const
 
G4double GetMeanExcitationEnergy () const
 
G4doubleGetShellCorrectionVector () const
 
G4double GetTau0 () const
 
G4double GetTaul () const
 
G4double GetZ () const
 
G4double GetZ3 () const
 
G4double GetZZ3 () const
 
G4bool operator!= (const G4IonisParamElm &) const =delete
 
G4IonisParamElmoperator= (const G4IonisParamElm &)=delete
 
G4bool operator== (const G4IonisParamElm &) const =delete
 
 ~G4IonisParamElm ()
 

Private Attributes

G4double fAlow
 
G4double fBetheBlochLow
 
G4double fBlow
 
G4double fClow
 
G4double fLFactor
 
G4double flogZ3
 
G4double fMeanExcitationEnergy
 
G4doublefShellCorrectionVector
 
G4double fTau0
 
G4double fTaul
 
G4double fVFermi
 
G4double fZ
 
G4double fZ3
 
G4double fZZ3
 

Detailed Description

Definition at line 50 of file G4IonisParamElm.hh.

Constructor & Destructor Documentation

◆ G4IonisParamElm() [1/3]

G4IonisParamElm::G4IonisParamElm ( G4double  Z)

Definition at line 51 of file G4IonisParamElm.cc.

52{
53 G4int Z = G4lrint(AtomNumber);
54 if (Z < 1) {
55 G4Exception("G4IonisParamElm::G4IonisParamElm()", "mat501", FatalException,
56 "It is not allowed to create an Element with Z<1");
57 }
58 G4Pow* g4pow = G4Pow::GetInstance();
59
60 // some basic functions of the atomic number
61 fZ = Z;
62 fZ3 = g4pow->Z13(Z);
63 fZZ3 = fZ3*g4pow->Z13(Z+1);
64 flogZ3 = g4pow->logZ(Z)/3.;
65
68
69 // compute parameters for ion transport
70 // The aproximation from:
71 // J.F.Ziegler, J.P. Biersack, U. Littmark
72 // The Stopping and Range of Ions in Matter,
73 // Vol.1, Pergamon Press, 1985
74 // Fast ions or hadrons
75
76 G4int iz = Z - 1;
77 if(91 < iz) { iz = 91; }
78
79 static const G4double vFermi[92] = {
80 1.0309, 0.15976, 0.59782, 1.0781, 1.0486, 1.0, 1.058, 0.93942, 0.74562, 0.3424,
81 0.45259, 0.71074, 0.90519, 0.97411, 0.97184, 0.89852, 0.70827, 0.39816, 0.36552, 0.62712,
82 0.81707, 0.9943, 1.1423, 1.2381, 1.1222, 0.92705, 1.0047, 1.2, 1.0661, 0.97411,
83 0.84912, 0.95, 1.0903, 1.0429, 0.49715, 0.37755, 0.35211, 0.57801, 0.77773, 1.0207,
84 1.029, 1.2542, 1.122, 1.1241, 1.0882, 1.2709, 1.2542, 0.90094, 0.74093, 0.86054,
85 0.93155, 1.0047, 0.55379, 0.43289, 0.32636, 0.5131, 0.695, 0.72591, 0.71202, 0.67413,
86 0.71418, 0.71453, 0.5911, 0.70263, 0.68049, 0.68203, 0.68121, 0.68532, 0.68715, 0.61884,
87 0.71801, 0.83048, 1.1222, 1.2381, 1.045, 1.0733, 1.0953, 1.2381, 1.2879, 0.78654,
88 0.66401, 0.84912, 0.88433, 0.80746, 0.43357, 0.41923, 0.43638, 0.51464, 0.73087, 0.81065,
89 1.9578, 1.0257} ;
90
91 static const G4double lFactor[92] = {
92 1.0, 1.0, 1.1, 1.06, 1.01, 1.03, 1.04, 0.99, 0.95, 0.9,
93 0.82, 0.81, 0.83, 0.88, 1.0, 0.95, 0.97, 0.99, 0.98, 0.97,
94 0.98, 0.97, 0.96, 0.93, 0.91, 0.9, 0.88, 0.9, 0.9, 0.9,
95 0.9, 0.85, 0.9, 0.9, 0.91, 0.92, 0.9, 0.9, 0.9, 0.9,
96 0.9, 0.88, 0.9, 0.88, 0.88, 0.9, 0.9, 0.88, 0.9, 0.9,
97 0.9, 0.9, 0.96, 1.2, 0.9, 0.88, 0.88, 0.85, 0.9, 0.9,
98 0.92, 0.95, 0.99, 1.03, 1.05, 1.07, 1.08, 1.1, 1.08, 1.08,
99 1.08, 1.08, 1.09, 1.09, 1.1, 1.11, 1.12, 1.13, 1.14, 1.15,
100 1.17, 1.2, 1.18, 1.17, 1.17, 1.16, 1.16, 1.16, 1.16, 1.16,
101 1.16, 1.16} ;
102
103 fVFermi = vFermi[iz];
104 fLFactor = lFactor[iz];
105
106 // obsolete parameters for ionisation
109
110 // compute the Bethe-Bloch formula for energy = fTaul*particle mass
112 G4double w = fTaul*(fTaul+2.) ;
113 fBetheBlochLow = (fTaul+1.)*(fTaul+1.)*std::log(2.*w/rate)/w - 1. ;
115
116 fClow = std::sqrt(fTaul)*fBetheBlochLow;
117 fAlow = 6.458040 * fClow/fTau0;
118 G4double Taum = 0.035*fZ3*MeV/proton_mass_c2;
119 fBlow =-3.229020*fClow/(fTau0*std::sqrt(Taum));
120
121 // Shell correction parameterization
122 fShellCorrectionVector = new G4double[3]; //[3]
123 rate = 0.001*fMeanExcitationEnergy/eV;
124 G4double rate2 = rate*rate;
125 /*
126 fShellCorrectionVector[0] = ( 1.10289e5 + 5.14781e8*rate)*rate2 ;
127 fShellCorrectionVector[1] = ( 7.93805e3 - 2.22565e7*rate)*rate2 ;
128 fShellCorrectionVector[2] = (-9.92256e1 + 2.10823e5*rate)*rate2 ;
129 */
130 fShellCorrectionVector[0] = ( 0.422377 + 3.858019*rate)*rate2 ;
131 fShellCorrectionVector[1] = ( 0.0304043 - 0.1667989*rate)*rate2 ;
132 fShellCorrectionVector[2] = (-0.00038106 + 0.00157955*rate)*rate2 ;
133}
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35
static constexpr double eV
Definition: G4SIunits.hh:201
static constexpr double MeV
Definition: G4SIunits.hh:200
double G4double
Definition: G4Types.hh:83
int G4int
Definition: G4Types.hh:85
const G4int Z[17]
G4double * fShellCorrectionVector
G4double fMeanExcitationEnergy
G4double GetMeanIonisationEnergy(G4int Z) const
static G4NistManager * Instance()
Definition: G4Pow.hh:49
static G4Pow * GetInstance()
Definition: G4Pow.cc:41
G4double logZ(G4int Z) const
Definition: G4Pow.hh:137
G4double Z13(G4int Z) const
Definition: G4Pow.hh:123
float electron_mass_c2
Definition: hepunit.py:273
int twopi_mc2_rcl2
Definition: hepunit.py:293
float proton_mass_c2
Definition: hepunit.py:274
int G4lrint(double ad)
Definition: templates.hh:134

References source.hepunit::electron_mass_c2, eV, fAlow, FatalException, fBetheBlochLow, fBlow, fClow, fLFactor, flogZ3, fMeanExcitationEnergy, fShellCorrectionVector, fTau0, fTaul, fVFermi, fZ, fZ3, fZZ3, G4Exception(), G4lrint(), G4Pow::GetInstance(), G4NistManager::GetMeanIonisationEnergy(), G4NistManager::Instance(), G4Pow::logZ(), MeV, source.hepunit::proton_mass_c2, source.hepunit::twopi_mc2_rcl2, Z, and G4Pow::Z13().

◆ ~G4IonisParamElm()

G4IonisParamElm::~G4IonisParamElm ( )

Definition at line 149 of file G4IonisParamElm.cc.

150{
152}

References fShellCorrectionVector.

◆ G4IonisParamElm() [2/3]

G4IonisParamElm::G4IonisParamElm ( __void__ &  )

◆ G4IonisParamElm() [3/3]

G4IonisParamElm::G4IonisParamElm ( G4IonisParamElm )
delete

Member Function Documentation

◆ GetAlow()

G4double G4IonisParamElm::GetAlow ( ) const
inline

Definition at line 69 of file G4IonisParamElm.hh.

69{return fAlow;}

References fAlow.

◆ GetBlow()

G4double G4IonisParamElm::GetBlow ( ) const
inline

Definition at line 70 of file G4IonisParamElm.hh.

70{return fBlow;}

References fBlow.

◆ GetClow()

G4double G4IonisParamElm::GetClow ( ) const
inline

Definition at line 71 of file G4IonisParamElm.hh.

71{return fClow;}

References fClow.

◆ GetFermiVelocity()

G4double G4IonisParamElm::GetFermiVelocity ( ) const
inline

Definition at line 77 of file G4IonisParamElm.hh.

77{return fVFermi;};

References fVFermi.

Referenced by G4IonisParamMat::ComputeIonParameters().

◆ GetLFactor()

G4double G4IonisParamElm::GetLFactor ( ) const
inline

Definition at line 78 of file G4IonisParamElm.hh.

78{return fLFactor;};

References fLFactor.

Referenced by G4IonisParamMat::ComputeIonParameters().

◆ GetlogZ3()

G4double G4IonisParamElm::GetlogZ3 ( ) const
inline

◆ GetMeanExcitationEnergy()

G4double G4IonisParamElm::GetMeanExcitationEnergy ( ) const
inline

Definition at line 74 of file G4IonisParamElm.hh.

References fMeanExcitationEnergy.

Referenced by G4IonisParamMat::ComputeMeanParameters().

◆ GetShellCorrectionVector()

G4double * G4IonisParamElm::GetShellCorrectionVector ( ) const
inline

Definition at line 80 of file G4IonisParamElm.hh.

References fShellCorrectionVector.

◆ GetTau0()

G4double G4IonisParamElm::GetTau0 ( ) const
inline

Definition at line 64 of file G4IonisParamElm.hh.

64{return fTau0;};

References fTau0.

◆ GetTaul()

G4double G4IonisParamElm::GetTaul ( ) const
inline

Definition at line 67 of file G4IonisParamElm.hh.

67{return fTaul;} // 2*MeV/proton mass

References fTaul.

◆ GetZ()

G4double G4IonisParamElm::GetZ ( ) const
inline

Definition at line 59 of file G4IonisParamElm.hh.

59{return fZ;} // atomic number Z

References fZ.

◆ GetZ3()

G4double G4IonisParamElm::GetZ3 ( ) const
inline

◆ GetZZ3()

G4double G4IonisParamElm::GetZZ3 ( ) const
inline

Definition at line 61 of file G4IonisParamElm.hh.

61{return fZZ3;} // std::pow (Z(Z+1),1/3)

References fZZ3.

◆ operator!=()

G4bool G4IonisParamElm::operator!= ( const G4IonisParamElm ) const
delete

◆ operator=()

G4IonisParamElm & G4IonisParamElm::operator= ( const G4IonisParamElm )
delete

◆ operator==()

G4bool G4IonisParamElm::operator== ( const G4IonisParamElm ) const
delete

Field Documentation

◆ fAlow

G4double G4IonisParamElm::fAlow
private

Definition at line 109 of file G4IonisParamElm.hh.

Referenced by G4IonisParamElm(), and GetAlow().

◆ fBetheBlochLow

G4double G4IonisParamElm::fBetheBlochLow
private

Definition at line 108 of file G4IonisParamElm.hh.

Referenced by G4IonisParamElm().

◆ fBlow

G4double G4IonisParamElm::fBlow
private

Definition at line 109 of file G4IonisParamElm.hh.

Referenced by G4IonisParamElm(), and GetBlow().

◆ fClow

G4double G4IonisParamElm::fClow
private

Definition at line 109 of file G4IonisParamElm.hh.

Referenced by G4IonisParamElm(), and GetClow().

◆ fLFactor

G4double G4IonisParamElm::fLFactor
private

Definition at line 115 of file G4IonisParamElm.hh.

Referenced by G4IonisParamElm(), and GetLFactor().

◆ flogZ3

G4double G4IonisParamElm::flogZ3
private

Definition at line 102 of file G4IonisParamElm.hh.

Referenced by G4IonisParamElm(), and GetlogZ3().

◆ fMeanExcitationEnergy

G4double G4IonisParamElm::fMeanExcitationEnergy
private

Definition at line 110 of file G4IonisParamElm.hh.

Referenced by G4IonisParamElm(), and GetMeanExcitationEnergy().

◆ fShellCorrectionVector

G4double* G4IonisParamElm::fShellCorrectionVector
private

Definition at line 111 of file G4IonisParamElm.hh.

Referenced by G4IonisParamElm(), GetShellCorrectionVector(), and ~G4IonisParamElm().

◆ fTau0

G4double G4IonisParamElm::fTau0
private

Definition at line 106 of file G4IonisParamElm.hh.

Referenced by G4IonisParamElm(), and GetTau0().

◆ fTaul

G4double G4IonisParamElm::fTaul
private

Definition at line 107 of file G4IonisParamElm.hh.

Referenced by G4IonisParamElm(), and GetTaul().

◆ fVFermi

G4double G4IonisParamElm::fVFermi
private

Definition at line 114 of file G4IonisParamElm.hh.

Referenced by G4IonisParamElm(), and GetFermiVelocity().

◆ fZ

G4double G4IonisParamElm::fZ
private

Definition at line 99 of file G4IonisParamElm.hh.

Referenced by G4IonisParamElm(), and GetZ().

◆ fZ3

G4double G4IonisParamElm::fZ3
private

Definition at line 100 of file G4IonisParamElm.hh.

Referenced by G4IonisParamElm(), and GetZ3().

◆ fZZ3

G4double G4IonisParamElm::fZZ3
private

Definition at line 101 of file G4IonisParamElm.hh.

Referenced by G4IonisParamElm(), and GetZZ3().


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