Geant4 11.3.0
Toolkit for the simulation of the passage of particles through matter
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Pages
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
 
std::size_t GetNumberOfIsotopes () const
 
G4IsotopeVectorGetIsotopeVector () const
 
G4doubleGetRelativeAbundanceVector () const
 
const G4IsotopeGetIsotope (G4int iso) const
 
std::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 const G4ElementTableGetElementTable ()
 
static std::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 66 of file G4Element.cc.

67 : fName(name), fSymbol(symbol)
68{
69 G4int iz = G4lrint(zeff);
70 if (iz < 1) {
72 ed << "Failed to create G4Element " << name << " Z= " << zeff << " < 1 !";
73 G4Exception("G4Element::G4Element()", "mat011", FatalException, ed);
74 }
75 if (std::abs(zeff - iz) > perMillion) {
77 ed << "G4Element Warning: " << name << " Z= " << zeff << " A= " << aeff / (g / mole);
78 G4Exception("G4Element::G4Element()", "mat017", JustWarning, ed);
79 }
80
81 InitializePointers();
82
83 fZeff = zeff;
84 fAeff = aeff;
85 fNeff = fAeff / (g / mole);
86
87 if (fNeff < 1.0) {
88 fNeff = 1.0;
89 }
90
91 if (fNeff < zeff) {
93 ed << "Failed to create G4Element " << name << " with Z= " << zeff << " N= " << fNeff
94 << " N < Z is not allowed" << G4endl;
95 G4Exception("G4Element::G4Element()", "mat012", FatalException, ed);
96 }
97
98 fNbOfAtomicShells = G4AtomicShells::GetNumberOfShells(iz);
99 fAtomicShells = new G4double[fNbOfAtomicShells];
100 fNbOfShellElectrons = new G4int[fNbOfAtomicShells];
101
102 AddNaturalIsotopes();
103
104 for (G4int i = 0; i < fNbOfAtomicShells; ++i) {
105 fAtomicShells[i] = G4AtomicShells::GetBindingEnergy(iz, i);
106 fNbOfShellElectrons[i] = G4AtomicShells::GetNumberOfElectrons(iz, i);
107 }
108 ComputeDerivedQuantities();
109}
@ 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

Referenced by G4Element(), GetElement(), operator!=(), operator<<, operator<<, operator=(), and operator==().

◆ G4Element() [2/3]

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

Definition at line 116 of file G4Element.cc.

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

◆ ~G4Element()

G4Element::~G4Element ( )
virtual

Definition at line 231 of file G4Element.cc.

232{
233 delete theIsotopeVector;
234 delete[] fRelativeAbundanceVector;
235 delete[] fAtomicShells;
236 delete[] fNbOfShellElectrons;
237 delete fIonisation;
238
239 // remove this element from the ElementTable
240 GetElementTableRef()[fIndexInTable] = nullptr;
241}

◆ G4Element() [3/3]

G4Element::G4Element ( G4Element & )
delete

Member Function Documentation

◆ AddIsotope()

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

Definition at line 138 of file G4Element.cc.

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

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

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

◆ GetElement()

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

Definition at line 412 of file G4Element.cc.

413{
414 // search the element by its name
415 for (auto const & J : GetElementTableRef())
416 {
417 if(J->GetName() == elementName)
418 {
419 return J;
420 }
421 }
422
423 // the element does not exist in the table
424 if (warning) {
425 G4cout << "\n---> warning from G4Element::GetElement(). The element: " << elementName
426 << " does not exist in the table. Return NULL pointer." << G4endl;
427 }
428 return nullptr;
429}
G4GLOB_DLL std::ostream G4cout

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

◆ GetElementTable()

const G4ElementTable * G4Element::GetElementTable ( )
static

Definition at line 401 of file G4Element.cc.

402{
403 return &GetElementTableRef();
404}

Referenced by G4FissLib::ApplyYourself(), G4NeutronHPCapture::ApplyYourself(), G4NeutronHPInelasticVI::ApplyYourself(), G4ParticleHPCaptureURR::ApplyYourself(), G4ParticleHPElastic::ApplyYourself(), G4ParticleHPElasticURR::ApplyYourself(), G4ParticleHPFission::ApplyYourself(), G4ParticleHPFissionURR::ApplyYourself(), G4ParticleHPInelastic::ApplyYourself(), G4ParticleHPInelasticURR::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(), 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(), G4ParticleHPProbabilityTablesStore::InitURRlimits(), G4NistManager::PrintG4Element(), and G4GDMLRead::StripNames().

◆ 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 223 of file G4Element.hh.

223{ return fNaturalAbundance; }

Referenced by G4CrossSectionDataStore::GetCrossSection().

◆ GetNbOfAtomicShells()

G4int G4Element::GetNbOfAtomicShells ( ) const
inline

◆ GetNbOfShellElectrons()

G4int G4Element::GetNbOfShellElectrons ( G4int index) const

Definition at line 372 of file G4Element.cc.

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

Referenced by G4KleinNishinaModel::SampleSecondaries().

◆ GetNumberOfElements()

◆ GetNumberOfIsotopes()

std::size_t G4Element::GetNumberOfIsotopes ( ) const
inline

◆ 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(), G4JAEAElasticScatteringModel::SampleSecondaries(), G4JAEAPolarizedElasticScatteringModel::SampleSecondaries(), G4LivermorePolarizedComptonModel::SampleSecondaries(), G4LowEPComptonModel::SampleSecondaries(), G4LowEPPolarizedComptonModel::SampleSecondaries(), G4MuBremsstrahlungModel::SampleSecondaries(), G4MuonToMuonPairProductionModel::SampleSecondaries(), G4MuPairProductionModel::SampleSecondaries(), G4PEEffectFluoModel::SampleSecondaries(), G4RiGeMuPairProductionModel::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 G4ParticleHPCaptureURR::ApplyYourself(), G4ParticleHPElasticURR::ApplyYourself(), G4ParticleHPFissionURR::ApplyYourself(), G4ParticleHPInelasticURR::ApplyYourself(), G4VCrossSectionDataSet::ComputeCrossSection(), G4NeutronCaptureXS::ComputeCrossSectionPerElement(), G4NeutronElasticXS::ComputeCrossSectionPerElement(), G4NeutronInelasticXS::ComputeCrossSectionPerElement(), G4ParticleInelasticXS::ComputeCrossSectionPerElement(), G4VCrossSectionDataSet::ComputeCrossSectionPerElement(), G4EmCalculator::ComputeCrossSectionPerShell(), G4CrossSectionDataStore::GetCrossSection(), G4NeutronHPCaptureData::GetCrossSection(), G4GammaConversionToMuons::GetCrossSectionPerAtom(), G4VComponentCrossSection::GetElasticElementCrossSection(), G4HadronicProcess::GetElementCrossSection(), G4ComponentGGHadronNucleusXsc::GetHadronNucleonXsc(), G4ComponentGGHadronNucleusXsc::GetHadronNucleonXscNS(), G4ComponentGGHadronNucleusXsc::GetHadronNucleonXscPDG(), G4ComponentGGHadronNucleusXsc::GetHNinelasticXsc(), G4VComponentCrossSection::GetInelasticElementCrossSection(), G4VComponentCrossSection::GetTotalElementCrossSection(), G4ParticleHPElementData::Init(), G4VEmModel::InitialiseForMaterial(), G4ElNeutrinoNucleusProcess::PostStepDoIt(), G4GammaConversionToMuons::PostStepDoIt(), G4MuNeutrinoNucleusProcess::PostStepDoIt(), G4TauNeutrinoNucleusProcess::PostStepDoIt(), G4AdjointBremsstrahlungModel::RapidSampleSecondaries(), G4BetheHeitler5DModel::SampleSecondaries(), G4BetheHeitlerModel::SampleSecondaries(), G4eBremsstrahlungRelModel::SampleSecondaries(), G4eCoulombScatteringModel::SampleSecondaries(), G4eDPWACoulombScatteringModel::SampleSecondaries(), G4hCoulombScatteringModel::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 225 of file G4Element.hh.

225{ fNaturalAbundance = val; }

Friends And Related Symbol Documentation

◆ operator<< [1/4]

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

Definition at line 458 of file G4Element.cc.

459{
460 flux << &element;
461 return flux;
462}

◆ operator<< [2/4]

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

Definition at line 433 of file G4Element.cc.

434{
435 std::ios::fmtflags mode = flux.flags();
436 flux.setf(std::ios::fixed, std::ios::floatfield);
437 G4long prec = flux.precision(3);
438
439 flux << " Element: " << element->fName << " (" << element->fSymbol << ")"
440 << " Z = " << std::setw(4) << std::setprecision(1) << element->fZeff
441 << " N = " << std::setw(5) << std::setprecision(1) << G4lrint(element->fNeff)
442 << " A = " << std::setw(6) << std::setprecision(3) << (element->fAeff) / (g / mole)
443 << " g/mole";
444
445 for (G4int i = 0; i < element->fNumberOfIsotopes; i++) {
446 flux << "\n ---> " << (*(element->theIsotopeVector))[i]
447 << " abundance: " << std::setw(6) << std::setprecision(3)
448 << (element->fRelativeAbundanceVector[i]) / perCent << " %";
449 }
450
451 flux.precision(prec);
452 flux.setf(mode, std::ios::floatfield);
453 return flux;
454}
long G4long
Definition G4Types.hh:87

◆ operator<< [3/4]

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

Definition at line 466 of file G4Element.cc.

467{
468 // Dump info for all known elements
469 flux << "\n***** Table : Nb of elements = " << ElementTable.size() << " *****\n" << G4endl;
470
471 for (auto const & i : ElementTable) {
472 flux << i << G4endl << G4endl;
473 }
474
475 return flux;
476}

◆ operator<< [4/4]

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

Definition at line 480 of file G4Element.cc.

481{
482 // Dump info for all known elements
483 flux << "\n***** Vector : Nb of elements = " << ElementVector.size() << " *****\n" << G4endl;
484
485 for (auto const & i : ElementVector) {
486 flux << i << G4endl << G4endl;
487 }
488
489 return flux;
490}

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