Geant4 11.1.1
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)
 
void AddIsotope (G4Isotope *isotope, G4double RelativeAbundance)
 
virtual ~G4Element ()
 
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
 
 G4Element (__void__ &)
 
void SetName (const G4String &name)
 
 G4Element (G4Element &)=delete
 
const G4Elementoperator= (const G4Element &)=delete
 
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 &, const G4Element *)
 
std::ostream & operator<< (std::ostream &, const G4Element &)
 
std::ostream & operator<< (std::ostream &, const G4ElementTable &)
 
std::ostream & operator<< (std::ostream &, const G4ElementVector &)
 

Detailed Description

Definition at line 97 of file G4Element.hh.

Constructor & Destructor Documentation

◆ G4Element() [1/4]

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

Definition at line 74 of file G4Element.cc.

76 : fName(name), fSymbol(symbol)
77{
78 G4int iz = G4lrint(zeff);
79 if (iz < 1) {
81 ed << "Failed to create G4Element " << name
82 << " Z= " << zeff << " < 1 !";
83 G4Exception ("G4Element::G4Element()", "mat011", FatalException, ed);
84 }
85 if (std::abs(zeff - iz) > perMillion) {
87 ed << "G4Element Warning: " << name << " Z= " << zeff
88 << " A= " << aeff/(g/mole);
89 G4Exception("G4Element::G4Element()", "mat017", JustWarning, ed);
90 }
91
92 InitializePointers();
93
94 fZeff = zeff;
95 fAeff = aeff;
96 fNeff = fAeff/(g/mole);
97
98 if(fNeff < 1.0)
99 {
100 fNeff = 1.0;
101 }
102
103 if (fNeff < zeff) {
105 ed << "Failed to create G4Element " << name
106 << " with Z= " << zeff << " N= " << fNeff
107 << " N < Z is not allowed" << G4endl;
108 G4Exception("G4Element::G4Element()", "mat012", FatalException, ed);
109 }
110
111 fNbOfAtomicShells = G4AtomicShells::GetNumberOfShells(iz);
112 fAtomicShells = new G4double[fNbOfAtomicShells];
113 fNbOfShellElectrons = new G4int[fNbOfAtomicShells];
114
115 AddNaturalIsotopes();
116
117 for (G4int i=0; i<fNbOfAtomicShells; ++i)
118 {
119 fAtomicShells[i] = G4AtomicShells::GetBindingEnergy(iz, i);
120 fNbOfShellElectrons[i] = G4AtomicShells::GetNumberOfElectrons(iz, i);
121 }
122 ComputeDerivedQuantities();
123}
@ JustWarning
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:59
std::ostringstream G4ExceptionDescription
Definition: G4Exception.hh:40
double G4double
Definition: G4Types.hh:83
int G4int
Definition: G4Types.hh:85
#define G4endl
Definition: G4ios.hh:57
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/4]

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

Definition at line 130 of file G4Element.cc.

132 : fName(name),fSymbol(symbol)
133{
134 InitializePointers();
135
136 size_t n = size_t(nIsotopes);
137
138 if(0 >= nIsotopes) {
140 ed << "Failed to create G4Element " << name
141 << " <" << symbol << "> with " << nIsotopes
142 << " isotopes.";
143 G4Exception ("G4Element::G4Element()", "mat012", FatalException, ed);
144 } else {
145 theIsotopeVector = new G4IsotopeVector(n, nullptr);
146 fRelativeAbundanceVector = new G4double[nIsotopes];
147 }
148}
std::vector< G4Isotope * > G4IsotopeVector

◆ ~G4Element()

G4Element::~G4Element ( )
virtual

Definition at line 256 of file G4Element.cc.

257{
258 delete theIsotopeVector;
259 delete[] fRelativeAbundanceVector;
260 delete[] fAtomicShells;
261 delete[] fNbOfShellElectrons;
262 delete fIonisation;
263
264 //remove this element from theElementTable
265 theElementTable[fIndexInTable] = nullptr;
266}

◆ G4Element() [3/4]

G4Element::G4Element ( __void__ &  )

Definition at line 248 of file G4Element.cc.

249 : fZeff(0), fNeff(0), fAeff(0)
250{
251 InitializePointers();
252}

◆ G4Element() [4/4]

G4Element::G4Element ( G4Element )
delete

Member Function Documentation

◆ AddIsotope()

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

Definition at line 154 of file G4Element.cc.

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

Referenced by G4GDMLReadMaterials::MixtureRead().

◆ GetA()

G4double G4Element::GetA ( ) const
inline

◆ GetAtomicMassAmu()

G4double G4Element::GetAtomicMassAmu ( ) const
inline

Definition at line 136 of file G4Element.hh.

136{return fNeff;}

◆ GetAtomicShell()

G4double G4Element::GetAtomicShell ( G4int  index) const

Definition at line 372 of file G4Element.cc.

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

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

◆ GetElement()

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

Definition at line 417 of file G4Element.cc.

418{
419 // search the element by its name
420 for(auto J : theElementTable)
421 {
422 if(J->GetName() == elementName)
423 {
424 return J;
425 }
426 }
427
428 // the element does not exist in the table
429 if (warning) {
430 G4cout << "\n---> warning from G4Element::GetElement(). The element: "
431 << elementName << " does not exist in the table. Return NULL pointer."
432 << G4endl;
433 }
434 return nullptr;
435}
G4GLOB_DLL std::ostream G4cout

Referenced by G4GDMLReadMaterials::GetElement().

◆ GetElementTable()

G4ElementTable * G4Element::GetElementTable ( )
static

Definition at line 403 of file G4Element.cc.

404{
405 return &theElementTable;
406}

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

◆ GetfCoulomb()

◆ GetfRadTsai()

G4double G4Element::GetfRadTsai ( ) const
inline

Definition at line 194 of file G4Element.hh.

194{return fRadTsai;}

◆ GetIndex()

◆ GetIonisation()

◆ GetIsotope()

◆ GetIsotopeVector()

◆ GetN()

◆ GetName()

◆ GetNaturalAbundanceFlag()

G4bool G4Element::GetNaturalAbundanceFlag ( ) const
inline

Definition at line 262 of file G4Element.hh.

263{
264 return fNaturalAbundance;
265}

Referenced by G4CrossSectionDataStore::GetCrossSection().

◆ GetNbOfAtomicShells()

G4int G4Element::GetNbOfAtomicShells ( ) const
inline

◆ GetNbOfShellElectrons()

G4int G4Element::GetNbOfShellElectrons ( G4int  index) const

Definition at line 387 of file G4Element.cc.

388{
389 if (i<0 || i>=fNbOfAtomicShells) {
391 ed << "Invalid argument " << i << " for G4Element " << fName
392 << " with Z= " << fZeff
393 << " and Nshells= " << fNbOfAtomicShells;
394 G4Exception("G4Element::GetNbOfShellElectrons()", "mat016",
395 FatalException, ed);
396 return 0;
397 }
398 return fNbOfShellElectrons[i];
399}

Referenced by G4KleinNishinaModel::SampleSecondaries().

◆ GetNumberOfElements()

◆ GetNumberOfIsotopes()

◆ GetRelativeAbundanceVector()

◆ GetSymbol()

const G4String & G4Element::GetSymbol ( ) const
inline

Definition at line 128 of file G4Element.hh.

128{return fSymbol;}

◆ GetZ()

G4double G4Element::GetZ ( ) const
inline

Definition at line 131 of file G4Element.hh.

131{return fZeff;}

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

◆ GetZasInt()

G4int G4Element::GetZasInt ( ) const
inline

Definition at line 132 of file G4Element.hh.

132{return fZ;}

Referenced by G4VCrossSectionDataSet::ComputeCrossSection(), G4NeutronCaptureXS::ComputeCrossSectionPerElement(), G4NeutronElasticXS::ComputeCrossSectionPerElement(), G4NeutronInelasticXS::ComputeCrossSectionPerElement(), G4ParticleInelasticXS::ComputeCrossSectionPerElement(), G4VCrossSectionDataSet::ComputeCrossSectionPerElement(), G4EmCalculator::ComputeCrossSectionPerShell(), G4IonisParamMat::ComputeDensityEffectOnFly(), G4CrossSectionDataStore::GetCrossSection(), G4GammaConversionToMuons::GetCrossSectionPerAtom(), G4VComponentCrossSection::GetElasticElementCrossSection(), G4HadronicProcess::GetElementCrossSection(), G4ComponentGGHadronNucleusXsc::GetHadronNucleonXsc(), G4ComponentGGHadronNucleusXsc::GetHadronNucleonXscNS(), G4ComponentGGHadronNucleusXsc::GetHadronNucleonXscPDG(), G4ComponentGGHadronNucleusXsc::GetHNinelasticXsc(), G4VComponentCrossSection::GetInelasticElementCrossSection(), G4VComponentCrossSection::GetTotalElementCrossSection(), G4VEmModel::InitialiseForMaterial(), G4GammaConversionToMuons::PostStepDoIt(), G4AdjointBremsstrahlungModel::RapidSampleSecondaries(), G4LivermoreBremsstrahlungModel::SampleSecondaries(), G4eBremsstrahlungRelModel::SampleSecondaries(), G4SeltzerBergerModel::SampleSecondaries(), G4eSingleCoulombScatteringModel::SampleSecondaries(), G4IonCoulombScatteringModel::SampleSecondaries(), G4LivermoreComptonModel::SampleSecondaries(), G4LivermorePhotoElectricModel::SampleSecondaries(), G4LivermorePolarizedRayleighModel::SampleSecondaries(), G4PenelopePhotoElectricModel::SampleSecondaries(), G4BetheHeitlerModel::SampleSecondaries(), G4eCoulombScatteringModel::SampleSecondaries(), G4eDPWACoulombScatteringModel::SampleSecondaries(), G4hCoulombScatteringModel::SampleSecondaries(), G4KleinNishinaModel::SampleSecondaries(), G4PairProductionRelModel::SampleSecondaries(), G4BetheHeitler5DModel::SampleSecondaries(), G4CrossSectionDataStore::SampleZandA(), 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 214 of file G4Element.hh.

214{fName=name;}

Referenced by G4GDMLRead::StripNames().

◆ SetNaturalAbundanceFlag()

void G4Element::SetNaturalAbundanceFlag ( G4bool  val)
inline

Definition at line 267 of file G4Element.hh.

268{
269 fNaturalAbundance = val;
270}

Friends And Related Function Documentation

◆ operator<< [1/4]

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

Definition at line 467 of file G4Element.cc.

468{
469 flux << &element;
470 return flux;
471}

◆ operator<< [2/4]

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

Definition at line 439 of file G4Element.cc.

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

◆ operator<< [3/4]

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

Definition at line 475 of file G4Element.cc.

476{
477 //Dump info for all known elements
478 flux << "\n***** Table : Nb of elements = " << ElementTable.size()
479 << " *****\n" << G4endl;
480
481 for(auto i : ElementTable)
482 {
483 flux << i << G4endl << G4endl;
484 }
485
486 return flux;
487}

◆ operator<< [4/4]

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

Definition at line 491 of file G4Element.cc.

492{
493 //Dump info for all known elements
494 flux << "\n***** Vector : Nb of elements = " << ElementVector.size()
495 << " *****\n" << G4endl;
496
497 for(auto i : ElementVector)
498 {
499 flux << i << G4endl << G4endl;
500 }
501
502 return flux;
503}

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