Geant4 11.2.2
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
G4Element Class Reference

#include <G4Element.hh>

Public Member Functions

 G4Element (const G4String &name, const G4String &symbol, G4double Zeff, G4double Aeff)
 
 G4Element (const G4String &name, const G4String &symbol, G4int nbIsotopes)
 
virtual ~G4Element ()
 
 G4Element (G4Element &)=delete
 
const G4Elementoperator= (const G4Element &)=delete
 
void AddIsotope (G4Isotope *isotope, G4double RelativeAbundance)
 
const G4StringGetName () const
 
const G4StringGetSymbol () const
 
G4double GetZ () const
 
G4int GetZasInt () const
 
G4double GetN () const
 
G4double GetAtomicMassAmu () const
 
G4double GetA () const
 
G4bool GetNaturalAbundanceFlag () const
 
void SetNaturalAbundanceFlag (G4bool)
 
G4int GetNbOfAtomicShells () const
 
G4double GetAtomicShell (G4int index) const
 
G4int GetNbOfShellElectrons (G4int index) const
 
size_t GetNumberOfIsotopes () const
 
G4IsotopeVectorGetIsotopeVector () const
 
G4doubleGetRelativeAbundanceVector () const
 
const G4IsotopeGetIsotope (G4int iso) const
 
size_t GetIndex () const
 
G4double GetfCoulomb () const
 
G4double GetfRadTsai () const
 
G4IonisParamElmGetIonisation () const
 
void SetName (const G4String &name)
 
G4bool operator== (const G4Element &) const =delete
 
G4bool operator!= (const G4Element &) const =delete
 

Static Public Member Functions

static G4ElementTableGetElementTable ()
 
static size_t GetNumberOfElements ()
 
static G4ElementGetElement (const G4String &name, G4bool warning=true)
 

Friends

std::ostream & operator<< (std::ostream &flux, const G4Element *element)
 
std::ostream & operator<< (std::ostream &flux, const G4Element &element)
 
std::ostream & operator<< (std::ostream &flux, const G4ElementTable &ElementTable)
 
std::ostream & operator<< (std::ostream &flux, const G4ElementVector &ElementVector)
 

Detailed Description

Definition at line 90 of file G4Element.hh.

Constructor & Destructor Documentation

◆ G4Element() [1/3]

G4Element::G4Element ( const G4String & name,
const G4String & symbol,
G4double Zeff,
G4double Aeff )

Definition at line 68 of file G4Element.cc.

69 : fName(name), fSymbol(symbol)
70{
71 G4int iz = G4lrint(zeff);
72 if (iz < 1) {
74 ed << "Failed to create G4Element " << name << " Z= " << zeff << " < 1 !";
75 G4Exception("G4Element::G4Element()", "mat011", FatalException, ed);
76 }
77 if (std::abs(zeff - iz) > perMillion) {
79 ed << "G4Element Warning: " << name << " Z= " << zeff << " A= " << aeff / (g / mole);
80 G4Exception("G4Element::G4Element()", "mat017", JustWarning, ed);
81 }
82
83 InitializePointers();
84
85 fZeff = zeff;
86 fAeff = aeff;
87 fNeff = fAeff / (g / mole);
88
89 if (fNeff < 1.0) {
90 fNeff = 1.0;
91 }
92
93 if (fNeff < zeff) {
95 ed << "Failed to create G4Element " << name << " with Z= " << zeff << " N= " << fNeff
96 << " N < Z is not allowed" << G4endl;
97 G4Exception("G4Element::G4Element()", "mat012", FatalException, ed);
98 }
99
100 fNbOfAtomicShells = G4AtomicShells::GetNumberOfShells(iz);
101 fAtomicShells = new G4double[fNbOfAtomicShells];
102 fNbOfShellElectrons = new G4int[fNbOfAtomicShells];
103
104 AddNaturalIsotopes();
105
106 for (G4int i = 0; i < fNbOfAtomicShells; ++i) {
107 fAtomicShells[i] = G4AtomicShells::GetBindingEnergy(iz, i);
108 fNbOfShellElectrons[i] = G4AtomicShells::GetNumberOfElectrons(iz, i);
109 }
110 ComputeDerivedQuantities();
111}
@ JustWarning
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
std::ostringstream G4ExceptionDescription
double G4double
Definition G4Types.hh:83
int G4int
Definition G4Types.hh:85
#define G4endl
Definition G4ios.hh:67
static G4int GetNumberOfElectrons(G4int Z, G4int SubshellNb)
static G4double GetBindingEnergy(G4int Z, G4int SubshellNb)
static G4int GetNumberOfShells(G4int Z)
const char * name(G4int ptype)
int G4lrint(double ad)
Definition templates.hh:134

◆ G4Element() [2/3]

G4Element::G4Element ( const G4String & name,
const G4String & symbol,
G4int nbIsotopes )

Definition at line 118 of file G4Element.cc.

119 : fName(name), fSymbol(symbol)
120{
121 InitializePointers();
122
123 auto n = size_t(nIsotopes);
124
125 if (0 >= nIsotopes) {
127 ed << "Failed to create G4Element " << name << " <" << symbol << "> with " << nIsotopes
128 << " isotopes.";
129 G4Exception("G4Element::G4Element()", "mat012", FatalException, ed);
130 }
131 else {
132 theIsotopeVector = new G4IsotopeVector(n, nullptr);
133 fRelativeAbundanceVector = new G4double[nIsotopes];
134 }
135}
std::vector< G4Isotope * > G4IsotopeVector

◆ ~G4Element()

G4Element::~G4Element ( )
virtual

Definition at line 234 of file G4Element.cc.

235{
236 delete theIsotopeVector;
237 delete[] fRelativeAbundanceVector;
238 delete[] fAtomicShells;
239 delete[] fNbOfShellElectrons;
240 delete fIonisation;
241
242 // remove this element from theElementTable
243 theElementTable[fIndexInTable] = nullptr;
244}

◆ G4Element() [3/3]

G4Element::G4Element ( G4Element & )
delete

Member Function Documentation

◆ AddIsotope()

void G4Element::AddIsotope ( G4Isotope * isotope,
G4double RelativeAbundance )

Definition at line 141 of file G4Element.cc.

142{
143 if (theIsotopeVector == nullptr) {
145 ed << "Failed to add Isotope to G4Element " << fName << " with Z= " << fZeff
146 << " N= " << fNeff;
147 G4Exception("G4Element::AddIsotope()", "mat013", FatalException, ed);
148 return;
149 }
150 G4int iz = isotope->GetZ();
151
152 // filling ...
153 if (fNumberOfIsotopes < (G4int)theIsotopeVector->size()) {
154 // check same Z
155 if (fNumberOfIsotopes == 0) {
156 fZeff = G4double(iz);
157 }
158 else if (G4double(iz) != fZeff) {
160 ed << "Failed to add Isotope Z= " << iz << " to G4Element " << fName
161 << " with different Z= " << fZeff << fNeff;
162 G4Exception("G4Element::AddIsotope()", "mat014", FatalException, ed);
163 return;
164 }
165 // Z ok
166 fRelativeAbundanceVector[fNumberOfIsotopes] = abundance;
167 (*theIsotopeVector)[fNumberOfIsotopes] = isotope;
168 ++fNumberOfIsotopes;
169 }
170 else {
172 ed << "Failed to add Isotope Z= " << iz << " to G4Element " << fName
173 << " - more isotopes than declared.";
174 G4Exception("G4Element::AddIsotope()", "mat015", FatalException, ed);
175 return;
176 }
177
178 // filled.
179 if (fNumberOfIsotopes == (G4int)theIsotopeVector->size()) {
180 G4double wtSum = 0.0;
181 fAeff = 0.0;
182 for (G4int i = 0; i < fNumberOfIsotopes; ++i) {
183 fAeff += fRelativeAbundanceVector[i] * (*theIsotopeVector)[i]->GetA();
184 wtSum += fRelativeAbundanceVector[i];
185 }
186 if (wtSum > 0.0) {
187 fAeff /= wtSum;
188 }
189 fNeff = fAeff / (g / mole);
190
191 if (wtSum != 1.0) {
192 for (G4int i = 0; i < fNumberOfIsotopes; ++i) {
193 fRelativeAbundanceVector[i] /= wtSum;
194 }
195 }
196
197 fNbOfAtomicShells = G4AtomicShells::GetNumberOfShells(iz);
198 fAtomicShells = new G4double[fNbOfAtomicShells];
199 fNbOfShellElectrons = new G4int[fNbOfAtomicShells];
200
201 for (G4int j = 0; j < fNbOfAtomicShells; ++j) {
202 fAtomicShells[j] = G4AtomicShells::GetBindingEnergy(iz, j);
203 fNbOfShellElectrons[j] = G4AtomicShells::GetNumberOfElectrons(iz, j);
204 }
205 ComputeDerivedQuantities();
206 }
207}
G4int GetZ() const
Definition G4Isotope.hh:80

Referenced by G4GDMLReadMaterials::MixtureRead().

◆ GetA()

G4double G4Element::GetA ( ) const
inline

◆ GetAtomicMassAmu()

G4double G4Element::GetAtomicMassAmu ( ) const
inline

Definition at line 124 of file G4Element.hh.

124{ return fNeff; }

Referenced by G4ChannelingFastSimCrystalData::SetMaterialProperties().

◆ GetAtomicShell()

G4double G4Element::GetAtomicShell ( G4int index) const

Definition at line 361 of file G4Element.cc.

362{
363 if (i < 0 || i >= fNbOfAtomicShells) {
365 ed << "Invalid argument " << i << " in for G4Element " << fName << " with Z= " << fZeff
366 << " and Nshells= " << fNbOfAtomicShells;
367 G4Exception("G4Element::GetAtomicShell()", "mat016", FatalException, ed);
368 return 0.0;
369 }
370 return fAtomicShells[i];
371}

Referenced by G4AdjointPhotoElectricModel::AdjointCrossSectionPerAtom(), G4KleinNishinaModel::SampleSecondaries(), and G4PEEffectFluoModel::SampleSecondaries().

◆ GetElement()

G4Element * G4Element::GetElement ( const G4String & name,
G4bool warning = true )
static

Definition at line 397 of file G4Element.cc.

398{
399 // search the element by its name
400 for (auto J : theElementTable) {
401 if (J->GetName() == elementName) {
402 return J;
403 }
404 }
405
406 // the element does not exist in the table
407 if (warning) {
408 G4cout << "\n---> warning from G4Element::GetElement(). The element: " << elementName
409 << " does not exist in the table. Return NULL pointer." << G4endl;
410 }
411 return nullptr;
412}
G4GLOB_DLL std::ostream G4cout

Referenced by G4GDMLReadMaterials::GetElement(), and G4XrayReflection::SaveHenkeDataAsMaterialProperty().

◆ GetElementTable()

G4ElementTable * G4Element::GetElementTable ( )
static

Definition at line 389 of file G4Element.cc.

389{ return &theElementTable; }

Referenced by G4FissLib::ApplyYourself(), G4NeutronHPCapture::ApplyYourself(), G4NeutronHPInelasticVI::ApplyYourself(), G4ParticleHPElastic::ApplyYourself(), G4ParticleHPFission::ApplyYourself(), G4ParticleHPInelastic::ApplyYourself(), G4AdjointCSManager::BuildCrossSectionMatrices(), G4KokoulinMuonNuclearXS::BuildCrossSectionTable(), G4CrossSectionHP::BuildPhysicsTable(), G4GammaNuclearXS::BuildPhysicsTable(), G4NeutronCaptureXS::BuildPhysicsTable(), G4NeutronElasticXS::BuildPhysicsTable(), G4NeutronHPCapture::BuildPhysicsTable(), G4NeutronHPCaptureData::BuildPhysicsTable(), G4NeutronInelasticXS::BuildPhysicsTable(), G4ParticleHPElastic::BuildPhysicsTable(), G4ParticleHPElasticData::BuildPhysicsTable(), G4ParticleHPFission::BuildPhysicsTable(), G4ParticleHPFissionData::BuildPhysicsTable(), G4ParticleHPInelastic::BuildPhysicsTable(), G4ParticleHPInelasticData::BuildPhysicsTable(), G4ParticleHPJENDLHEData::BuildPhysicsTable(), G4ParticleHPThermalScatteringData::BuildPhysicsTable(), G4ParticleInelasticXS::BuildPhysicsTable(), G4HadronHElasticPhysics::ConstructProcess(), G4LENDCrossSection::create_used_target_map(), G4LENDModel::create_used_target_map(), G4CrossSectionHP::DumpPhysicsTable(), G4NeutronHPCaptureData::DumpPhysicsTable(), G4ParticleHPElasticData::DumpPhysicsTable(), G4ParticleHPFissionData::DumpPhysicsTable(), G4ParticleHPInelasticData::DumpPhysicsTable(), G4NistElementBuilder::FindElement(), G4NistElementBuilder::FindOrBuildElement(), G4FissLib::G4FissLib(), G4ParticleHPData::G4ParticleHPData(), G4NistManager::GetElement(), G4DiffuseElastic::Initialise(), G4DiffuseElasticV2::Initialise(), G4LivermoreBremsstrahlungModel::Initialise(), G4LivermoreComptonModel::Initialise(), G4LivermoreGammaConversion5DModel::Initialise(), G4LivermoreGammaConversionModel::Initialise(), G4LivermorePhotoElectricModel::Initialise(), G4LivermorePolarizedRayleighModel::Initialise(), G4LivermoreRayleighModel::Initialise(), G4NuclNuclDiffuseElastic::Initialise(), G4SeltzerBergerModel::Initialise(), G4VAtomDeexcitation::InitialiseAtomicDeexcitation(), G4BetheHeitlerModel::InitialiseElementData(), G4NistManager::PrintG4Element(), G4GDMLRead::StripNames(), and G4NistManager::~G4NistManager().

◆ GetfCoulomb()

◆ GetfRadTsai()

G4double G4Element::GetfRadTsai ( ) const
inline

Definition at line 168 of file G4Element.hh.

168{ return fRadTsai; }

◆ GetIndex()

◆ GetIonisation()

◆ GetIsotope()

◆ GetIsotopeVector()

◆ GetN()

◆ GetName()

◆ GetNaturalAbundanceFlag()

G4bool G4Element::GetNaturalAbundanceFlag ( ) const
inline

Definition at line 221 of file G4Element.hh.

221{ return fNaturalAbundance; }

Referenced by G4CrossSectionDataStore::GetCrossSection().

◆ GetNbOfAtomicShells()

G4int G4Element::GetNbOfAtomicShells ( ) const
inline

◆ GetNbOfShellElectrons()

G4int G4Element::GetNbOfShellElectrons ( G4int index) const

Definition at line 375 of file G4Element.cc.

376{
377 if (i < 0 || i >= fNbOfAtomicShells) {
379 ed << "Invalid argument " << i << " for G4Element " << fName << " with Z= " << fZeff
380 << " and Nshells= " << fNbOfAtomicShells;
381 G4Exception("G4Element::GetNbOfShellElectrons()", "mat016", FatalException, ed);
382 return 0;
383 }
384 return fNbOfShellElectrons[i];
385}

Referenced by G4KleinNishinaModel::SampleSecondaries().

◆ GetNumberOfElements()

◆ GetNumberOfIsotopes()

◆ GetRelativeAbundanceVector()

◆ GetSymbol()

const G4String & G4Element::GetSymbol ( ) const
inline

Definition at line 116 of file G4Element.hh.

116{ return fSymbol; }

◆ GetZ()

G4double G4Element::GetZ ( ) const
inline

Definition at line 119 of file G4Element.hh.

119{ return fZeff; }

Referenced by G4AdjointPhotoElectricModel::AdjointCrossSectionPerAtom(), G4ParticleHPElastic::ApplyYourself(), G4ParticleHPThermalScattering::ApplyYourself(), G4Nucleus::ChooseParameters(), G4AdjointCSManager::ComputeAdjointCS(), G4EmCalculator::ComputeCrossSectionPerAtom(), G4VEmModel::ComputeCrossSectionPerAtom(), G4ICRU49NuclearStoppingModel::ComputeDEDXPerVolume(), G4PAIxSection::ComputeLowEnergyCof(), G4PAIySection::ComputeLowEnergyCof(), G4HadronHElasticPhysics::ConstructProcess(), G4LENDCrossSection::create_used_target_map(), G4LENDModel::create_used_target_map(), G4GDMLWriteMaterials::ElementWrite(), G4ParticleHPElasticData::GetCrossSection(), G4ParticleHPFissionData::GetCrossSection(), G4ParticleHPInelasticData::GetCrossSection(), G4ParticleHPJENDLHEData::GetCrossSection(), GVFlashShowerParameterisation::GetEffZ(), G4ChargeExchangeProcess::GetElementCrossSection(), G4ParticleHPThermalBoost::GetThermalEnergy(), G4EmUtility::SampleRandomElement(), G4eBremParametrizedModel::SampleSecondaries(), G4IonCoulombScatteringModel::SampleSecondaries(), G4JAEAElasticScatteringModel::SampleSecondaries(), G4JAEAPolarizedElasticScatteringModel::SampleSecondaries(), G4LivermorePolarizedComptonModel::SampleSecondaries(), G4LowEPComptonModel::SampleSecondaries(), G4LowEPPolarizedComptonModel::SampleSecondaries(), G4MuBremsstrahlungModel::SampleSecondaries(), G4MuonToMuonPairProductionModel::SampleSecondaries(), G4MuPairProductionModel::SampleSecondaries(), G4PEEffectFluoModel::SampleSecondaries(), G4XrayRayleighModel::SampleSecondaries(), and G4ChannelingFastSimCrystalData::SetMaterialProperties().

◆ GetZasInt()

G4int G4Element::GetZasInt ( ) const
inline

Definition at line 120 of file G4Element.hh.

120{ return fZ; }

Referenced by G4VCrossSectionDataSet::ComputeCrossSection(), G4NeutronCaptureXS::ComputeCrossSectionPerElement(), G4NeutronElasticXS::ComputeCrossSectionPerElement(), G4NeutronInelasticXS::ComputeCrossSectionPerElement(), G4ParticleInelasticXS::ComputeCrossSectionPerElement(), G4VCrossSectionDataSet::ComputeCrossSectionPerElement(), G4EmCalculator::ComputeCrossSectionPerShell(), G4IonisParamMat::ComputeDensityEffectOnFly(), G4DensityEffectCalculator::G4DensityEffectCalculator(), G4CrossSectionDataStore::GetCrossSection(), G4NeutronHPCaptureData::GetCrossSection(), G4GammaConversionToMuons::GetCrossSectionPerAtom(), G4VComponentCrossSection::GetElasticElementCrossSection(), G4HadronicProcess::GetElementCrossSection(), G4ComponentGGHadronNucleusXsc::GetHadronNucleonXsc(), G4ComponentGGHadronNucleusXsc::GetHadronNucleonXscNS(), G4ComponentGGHadronNucleusXsc::GetHadronNucleonXscPDG(), G4ComponentGGHadronNucleusXsc::GetHNinelasticXsc(), G4VComponentCrossSection::GetInelasticElementCrossSection(), G4VComponentCrossSection::GetTotalElementCrossSection(), G4ParticleHPElementData::Init(), G4LivermorePhotoElectricModel::Initialise(), G4VEmModel::InitialiseForMaterial(), G4ElNeutrinoNucleusProcess::PostStepDoIt(), G4GammaConversionToMuons::PostStepDoIt(), G4MuNeutrinoNucleusProcess::PostStepDoIt(), G4TauNeutrinoNucleusProcess::PostStepDoIt(), G4AdjointBremsstrahlungModel::RapidSampleSecondaries(), G4ParticleHPChannel::Register(), G4BetheHeitler5DModel::SampleSecondaries(), G4BetheHeitlerModel::SampleSecondaries(), G4eBremsstrahlungRelModel::SampleSecondaries(), G4eCoulombScatteringModel::SampleSecondaries(), G4eDPWACoulombScatteringModel::SampleSecondaries(), G4eSingleCoulombScatteringModel::SampleSecondaries(), G4hCoulombScatteringModel::SampleSecondaries(), G4IonCoulombScatteringModel::SampleSecondaries(), G4KleinNishinaModel::SampleSecondaries(), G4LivermoreBremsstrahlungModel::SampleSecondaries(), G4LivermoreComptonModel::SampleSecondaries(), G4LivermorePhotoElectricModel::SampleSecondaries(), G4LivermorePolarizedRayleighModel::SampleSecondaries(), G4LivermoreRayleighModel::SampleSecondaries(), G4PAIModel::SampleSecondaries(), G4PAIPhotModel::SampleSecondaries(), G4PairProductionRelModel::SampleSecondaries(), G4PenelopePhotoElectricModel::SampleSecondaries(), G4SeltzerBergerModel::SampleSecondaries(), G4CrossSectionDataStore::SampleZandA(), G4CrossSectionHP::SelectIsotope(), G4GammaNuclearXS::SelectIsotope(), G4NeutronCaptureXS::SelectIsotope(), G4NeutronInelasticXS::SelectIsotope(), G4ParticleInelasticXS::SelectIsotope(), G4VEmModel::SelectRandomAtomNumber(), and G4ElementSelector::SelectZandA().

◆ operator!=()

G4bool G4Element::operator!= ( const G4Element & ) const
delete

◆ operator=()

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

◆ operator==()

G4bool G4Element::operator== ( const G4Element & ) const
delete

◆ SetName()

void G4Element::SetName ( const G4String & name)
inline

Definition at line 179 of file G4Element.hh.

179{ fName = name; }

Referenced by G4GDMLRead::StripNames().

◆ SetNaturalAbundanceFlag()

void G4Element::SetNaturalAbundanceFlag ( G4bool val)
inline

Definition at line 223 of file G4Element.hh.

223{ fNaturalAbundance = val; }

Friends And Related Symbol Documentation

◆ operator<< [1/4]

std::ostream & operator<< ( std::ostream & flux,
const G4Element & element )
friend

Definition at line 441 of file G4Element.cc.

442{
443 flux << &element;
444 return flux;
445}

◆ operator<< [2/4]

std::ostream & operator<< ( std::ostream & flux,
const G4Element * element )
friend

Definition at line 416 of file G4Element.cc.

417{
418 std::ios::fmtflags mode = flux.flags();
419 flux.setf(std::ios::fixed, std::ios::floatfield);
420 G4long prec = flux.precision(3);
421
422 flux << " Element: " << element->fName << " (" << element->fSymbol << ")"
423 << " Z = " << std::setw(4) << std::setprecision(1) << element->fZeff
424 << " N = " << std::setw(5) << std::setprecision(1) << G4lrint(element->fNeff)
425 << " A = " << std::setw(6) << std::setprecision(3) << (element->fAeff) / (g / mole)
426 << " g/mole";
427
428 for (G4int i = 0; i < element->fNumberOfIsotopes; i++) {
429 flux << "\n ---> " << (*(element->theIsotopeVector))[i]
430 << " abundance: " << std::setw(6) << std::setprecision(3)
431 << (element->fRelativeAbundanceVector[i]) / perCent << " %";
432 }
433
434 flux.precision(prec);
435 flux.setf(mode, std::ios::floatfield);
436 return flux;
437}
long G4long
Definition G4Types.hh:87

◆ operator<< [3/4]

std::ostream & operator<< ( std::ostream & flux,
const G4ElementTable & ElementTable )
friend

Definition at line 449 of file G4Element.cc.

450{
451 // Dump info for all known elements
452 flux << "\n***** Table : Nb of elements = " << ElementTable.size() << " *****\n" << G4endl;
453
454 for (auto i : ElementTable) {
455 flux << i << G4endl << G4endl;
456 }
457
458 return flux;
459}

◆ operator<< [4/4]

std::ostream & operator<< ( std::ostream & flux,
const G4ElementVector & ElementVector )
friend

Definition at line 463 of file G4Element.cc.

464{
465 // Dump info for all known elements
466 flux << "\n***** Vector : Nb of elements = " << ElementVector.size() << " *****\n" << G4endl;
467
468 for (auto i : ElementVector) {
469 flux << i << G4endl << G4endl;
470 }
471
472 return flux;
473}

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