Geant4 10.7.0
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)
 

Static Public Member Functions

static G4ElementTableGetElementTable ()
 
static size_t GetNumberOfElements ()
 
static G4ElementGetElement (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 &, G4ElementTable)
 

Detailed Description

Definition at line 96 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 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) fNeff = 1.0;
99
100 if (fNeff < zeff) {
102 ed << "Failed to create G4Element " << name
103 << " with Z= " << zeff << " N= " << fNeff
104 << " N < Z is not allowed" << G4endl;
105 G4Exception("G4Element::G4Element()", "mat012", FatalException, ed);
106 }
107
108 fNbOfAtomicShells = G4AtomicShells::GetNumberOfShells(iz);
109 fAtomicShells = new G4double[fNbOfAtomicShells];
110 fNbOfShellElectrons = new G4int[fNbOfAtomicShells];
111
112 AddNaturalIsotopes();
113
114 for (G4int i=0;i<fNbOfAtomicShells;i++)
115 {
116 fAtomicShells[i] = G4AtomicShells::GetBindingEnergy(iz, i);
117 fNbOfShellElectrons[i] = G4AtomicShells::GetNumberOfElectrons(iz, i);
118 }
119 ComputeDerivedQuantities();
120}
@ JustWarning
@ FatalException
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35
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/3]

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

Definition at line 127 of file G4Element.cc.

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

◆ ~G4Element()

G4Element::~G4Element ( )
virtual

Definition at line 252 of file G4Element.cc.

253{
254 if (theIsotopeVector) { delete theIsotopeVector; }
255 if (fRelativeAbundanceVector) { delete [] fRelativeAbundanceVector; }
256 if (fAtomicShells) { delete [] fAtomicShells; }
257 if (fNbOfShellElectrons) { delete [] fNbOfShellElectrons; }
258 if (fIonisation) { delete fIonisation; }
259
260 //remove this element from theElementTable
261 theElementTable[fIndexInTable] = 0;
262}

◆ G4Element() [3/3]

G4Element::G4Element ( __void__ &  )

Definition at line 244 of file G4Element.cc.

245 : fZeff(0), fNeff(0), fAeff(0)
246{
247 InitializePointers();
248}

Member Function Documentation

◆ AddIsotope()

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

Definition at line 151 of file G4Element.cc.

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

Referenced by G4tgbElement::BuildG4ElementFromIsotopes(), and G4GDMLReadMaterials::MixtureRead().

◆ GetA()

G4double G4Element::GetA ( ) const
inline

◆ GetAtomicMassAmu()

G4double G4Element::GetAtomicMassAmu ( ) const
inline

Definition at line 135 of file G4Element.hh.

135{return fNeff;}

◆ GetAtomicShell()

G4double G4Element::GetAtomicShell ( G4int  index) const

Definition at line 366 of file G4Element.cc.

367{
368 if (i<0 || i>=fNbOfAtomicShells) {
370 ed << "Invalid argument " << i << " in for G4Element " << fName
371 << " with Z= " << fZeff
372 << " and Nshells= " << fNbOfAtomicShells;
373 G4Exception("G4Element::GetAtomicShell()", "mat016", FatalException, ed);
374 return 0.0;
375 }
376 return fAtomicShells[i];
377}

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

◆ GetElement()

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

Definition at line 411 of file G4Element.cc.

412{
413 // search the element by its name
414 for (size_t J=0; J<theElementTable.size(); ++J)
415 {
416 if (theElementTable[J]->GetName() == elementName)
417 return theElementTable[J];
418 }
419
420 // the element does not exist in the table
421 if (warning) {
422 G4cout << "\n---> warning from G4Element::GetElement(). The element: "
423 << elementName << " does not exist in the table. Return NULL pointer."
424 << G4endl;
425 }
426 return nullptr;
427}
G4GLOB_DLL std::ostream G4cout
const G4String & GetName() const
Definition: G4Element.hh:126

Referenced by G4GDMLReadMaterials::GetElement().

◆ GetElementTable()

G4ElementTable * G4Element::GetElementTable ( )
static

Definition at line 397 of file G4Element.cc.

398{
399 return &theElementTable;
400}

Referenced by G4FissLib::ApplyYourself(), G4ParticleHPCapture::ApplyYourself(), G4ParticleHPElastic::ApplyYourself(), G4ParticleHPFission::ApplyYourself(), G4ParticleHPInelastic::ApplyYourself(), G4AdjointCSManager::BuildCrossSectionMatrices(), G4KokoulinMuonNuclearXS::BuildCrossSectionTable(), G4VRangeToEnergyConverter::BuildLossTable(), G4ParticleHPCapture::BuildPhysicsTable(), G4ParticleHPCaptureData::BuildPhysicsTable(), G4ParticleHPElastic::BuildPhysicsTable(), G4ParticleHPElasticData::BuildPhysicsTable(), G4ParticleHPFission::BuildPhysicsTable(), G4ParticleHPFissionData::BuildPhysicsTable(), G4ParticleHPInelastic::BuildPhysicsTable(), G4ParticleHPInelasticData::BuildPhysicsTable(), G4ParticleHPJENDLHEData::BuildPhysicsTable(), G4ParticleHPThermalScatteringData::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(), 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()

G4IsotopeVector * G4Element::GetIsotopeVector ( ) const
inline

◆ GetN()

◆ GetName()

◆ GetNaturalAbundanceFlag()

G4bool G4Element::GetNaturalAbundanceFlag ( ) const
inline

Definition at line 261 of file G4Element.hh.

262{
263 return fNaturalAbundance;
264}

Referenced by G4CrossSectionDataStore::GetCrossSection().

◆ GetNbOfAtomicShells()

G4int G4Element::GetNbOfAtomicShells ( ) const
inline

◆ GetNbOfShellElectrons()

G4int G4Element::GetNbOfShellElectrons ( G4int  index) const

Definition at line 381 of file G4Element.cc.

382{
383 if (i<0 || i>=fNbOfAtomicShells) {
385 ed << "Invalid argument " << i << " for G4Element " << fName
386 << " with Z= " << fZeff
387 << " and Nshells= " << fNbOfAtomicShells;
388 G4Exception("G4Element::GetNbOfShellElectrons()", "mat016",
389 FatalException, ed);
390 return 0;
391 }
392 return fNbOfShellElectrons[i];
393}

Referenced by G4KleinNishinaModel::SampleSecondaries().

◆ GetNumberOfElements()

◆ GetNumberOfIsotopes()

◆ GetRelativeAbundanceVector()

◆ GetSymbol()

const G4String & G4Element::GetSymbol ( ) const
inline

Definition at line 127 of file G4Element.hh.

127{return fSymbol;}

◆ GetZ()

G4double G4Element::GetZ ( ) const
inline

Definition at line 130 of file G4Element.hh.

130{return fZeff;}

Referenced by G4VCrossSectionHandler::ActiveElements(), G4AdjointPhotoElectricModel::AdjointCrossSectionPerAtom(), G4ParticleHPThermalScattering::ApplyYourself(), G4AugerData::BuildAugerTransitionTable(), G4AdjointCSManager::BuildCrossSectionMatrices(), 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(), G4NeutrinoElectronProcess::PostStepDoIt(), G4ElNeutrinoNucleusProcess::PostStepDoIt(), G4MuNeutrinoNucleusProcess::PostStepDoIt(), G4ParticleHPChannel::Register(), G4eBremParametrizedModel::SampleSecondaries(), G4JAEAElasticScatteringModel::SampleSecondaries(), G4JAEAPolarizedElasticScatteringModel::SampleSecondaries(), G4LivermoreComptonModel::SampleSecondaries(), G4LivermoreComptonModifiedModel::SampleSecondaries(), G4LivermorePolarizedComptonModel::SampleSecondaries(), G4LivermorePolarizedPhotoElectricGDModel::SampleSecondaries(), G4LivermorePolarizedPhotoElectricModel::SampleSecondaries(), G4LivermorePolarizedRayleighModel::SampleSecondaries(), G4LivermoreRayleighModel::SampleSecondaries(), G4LowEPComptonModel::SampleSecondaries(), G4LowEPPolarizedComptonModel::SampleSecondaries(), G4XrayRayleighModel::SampleSecondaries(), G4IonCoulombScatteringModel::SampleSecondaries(), G4PAIModel::SampleSecondaries(), G4PAIPhotModel::SampleSecondaries(), G4MuBremsstrahlungModel::SampleSecondaries(), G4MuPairProductionModel::SampleSecondaries(), G4PEEffectFluoModel::SampleSecondaries(), and G4ElementSelector::SelectZandA().

◆ GetZasInt()

G4int G4Element::GetZasInt ( ) const
inline

Definition at line 131 of file G4Element.hh.

131{return fZ;}

Referenced by G4VCrossSectionDataSet::ComputeCrossSection(), G4EmCalculator::ComputeCrossSectionPerShell(), G4IonisParamMat::ComputeDensityEffectOnFly(), G4DensityEffectCalculator::G4DensityEffectCalculator(), G4CrossSectionDataStore::GetCrossSection(), G4GammaConversionToMuons::GetCrossSectionPerAtom(), G4VComponentCrossSection::GetElasticElementCrossSection(), G4HadronicProcess::GetElementCrossSection(), G4ComponentGGHadronNucleusXsc::GetHadronNucleonXsc(), G4ComponentGGHadronNucleusXsc::GetHadronNucleonXscNS(), G4ComponentGGHadronNucleusXsc::GetHadronNucleonXscPDG(), G4ComponentGGHadronNucleusXsc::GetHNinelasticXsc(), G4VComponentCrossSection::GetInelasticElementCrossSection(), G4VComponentCrossSection::GetTotalElementCrossSection(), G4BetheHeitlerModel::InitialiseElementData(), G4VEmModel::InitialiseForMaterial(), G4GammaConversionToMuons::PostStepDoIt(), G4AdjointBremsstrahlungModel::RapidSampleSecondaries(), G4LivermoreBremsstrahlungModel::SampleSecondaries(), G4eBremsstrahlungRelModel::SampleSecondaries(), G4SeltzerBergerModel::SampleSecondaries(), G4LivermorePhotoElectricModel::SampleSecondaries(), G4PenelopePhotoElectricModel::SampleSecondaries(), G4eSingleCoulombScatteringModel::SampleSecondaries(), G4IonCoulombScatteringModel::SampleSecondaries(), G4BetheHeitlerModel::SampleSecondaries(), G4eCoulombScatteringModel::SampleSecondaries(), G4eDPWACoulombScatteringModel::SampleSecondaries(), G4hCoulombScatteringModel::SampleSecondaries(), G4KleinNishinaModel::SampleSecondaries(), G4PairProductionRelModel::SampleSecondaries(), G4BetheHeitler5DModel::SampleSecondaries(), G4CrossSectionDataStore::SampleZandA(), G4NeutronCaptureXS::SelectIsotope(), G4NeutronInelasticXS::SelectIsotope(), G4ParticleInelasticXS::SelectIsotope(), and G4VEmModel::SelectRandomAtomNumber().

◆ SetName()

void G4Element::SetName ( const G4String name)
inline

Definition at line 213 of file G4Element.hh.

213{fName=name;}

Referenced by G4GDMLRead::StripNames().

◆ SetNaturalAbundanceFlag()

void G4Element::SetNaturalAbundanceFlag ( G4bool  val)
inline

Definition at line 266 of file G4Element.hh.

267{
268 fNaturalAbundance = val;
269}

Friends And Related Function Documentation

◆ operator<< [1/3]

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/3]

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

Definition at line 431 of file G4Element.cc.

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

◆ operator<< [3/3]

std::ostream & operator<< ( std::ostream &  flux,
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()
470 << " *****\n" << G4endl;
471
472 for (size_t i=0; i<ElementTable.size(); i++) flux << ElementTable[i]
473 << G4endl << G4endl;
474
475 return flux;
476}

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