Geant4  10.03.p03
 All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Groups Pages
G4Material.cc
Go to the documentation of this file.
1 //
2 // ********************************************************************
3 // * License and Disclaimer *
4 // * *
5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. *
10 // * *
11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitation of liability. *
17 // * *
18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Software license. *
24 // ********************************************************************
25 //
26 // $Id: G4Material.cc 100662 2016-10-31 10:21:14Z gcosmo $
27 //
28 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
29 //
30 // 26-06-96, Code uses operators (+=, *=, ++, -> etc.) correctly, P. Urban
31 // 10-07-96, new data members added by L.Urban
32 // 12-12-96, new data members added by L.Urban
33 // 20-01-97, aesthetic rearrangement. RadLength calculation modified.
34 // Data members Zeff and Aeff REMOVED (i.e. passed to the Elements).
35 // (local definition of Zeff in DensityEffect and FluctModel...)
36 // Vacuum defined as a G4State. Mixture flag removed, M.Maire.
37 // 29-01-97, State=Vacuum automatically set density=0 in the contructors.
38 // Subsequent protections have been put in the calculation of
39 // MeanExcEnergy, ShellCorrectionVector, DensityEffect, M.Maire.
40 // 11-02-97, ComputeDensityEffect() rearranged, M.Maire.
41 // 20-03-97, corrected initialization of pointers, M.Maire.
42 // 28-05-98, the kState=kVacuum has been removed.
43 // automatic check for a minimal density, M.Maire
44 // 12-06-98, new method AddMaterial() allowing mixture of materials, M.Maire
45 // 09-07-98, ionisation parameters removed from the class, M.Maire
46 // 05-10-98, change names: NumDensity -> NbOfAtomsPerVolume
47 // 18-11-98, new interface to SandiaTable
48 // 19-01-99 enlarge tolerance on test of coherence of gas conditions
49 // 19-07-99, Constructors with chemicalFormula added by V.Ivanchenko
50 // 16-01-01, Nuclear interaction length, M.Maire
51 // 12-03-01, G4bool fImplicitElement;
52 // copy constructor and assignement operator revised (mma)
53 // 03-05-01, flux.precision(prec) at begin/end of operator<<
54 // 17-07-01, migration to STL. M. Verderi.
55 // 14-09-01, Suppression of the data member fIndexInTable
56 // 26-02-02, fIndexInTable renewed
57 // 16-04-02, G4Exception put in constructor with chemical formula
58 // 06-05-02, remove the check of the ideal gas state equation
59 // 06-08-02, remove constructors with chemical formula (mma)
60 // 22-01-04, proper STL handling of theElementVector (Hisaya)
61 // 30-03-05, warning in GetMaterial(materialName)
62 // 09-03-06, minor change of printout (V.Ivanchenko)
63 // 10-01-07, compute fAtomVector in the case of mass fraction (V.Ivanchenko)
64 // 27-07-07, improve destructor (V.Ivanchenko)
65 // 18-10-07, move definition of material index to InitialisePointers (V.Ivanchenko)
66 // 13-08-08, do not use fixed size arrays (V.Ivanchenko)
67 // 26-10-11, new scheme for G4Exception (mma)
68 // 13-04-12, map<G4Material*,G4double> fMatComponents, filled in AddMaterial()
69 // 21-04-12, fMassOfMolecule, computed for AtomsCount (mma)
70 //
71 //
72 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
73 
74 #include <iomanip>
75 
76 #include "G4Material.hh"
77 #include "G4NistManager.hh"
78 #include "G4UnitsTable.hh"
79 #include "G4PhysicalConstants.hh"
80 #include "G4SystemOfUnits.hh"
81 #include "G4Exp.hh"
82 #include "G4Log.hh"
83 
84 G4MaterialTable G4Material::theMaterialTable;
85 
86 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
87 
88 // Constructor to create a material from scratch
89 
91  G4double a, G4double density,
92  G4State state, G4double temp, G4double pressure)
93  : fName(name)
94 {
95  InitializePointers();
96 
97  if (density < universe_mean_density)
98  {
99  G4cout << " G4Material WARNING:"
100  << " define a material with density=0 is not allowed. \n"
101  << " The material " << name << " will be constructed with the"
102  << " default minimal density: " << universe_mean_density/(g/cm3)
103  << "g/cm3" << G4endl;
104  density = universe_mean_density;
105  }
106 
107  fDensity = density;
108  fState = state;
109  fTemp = temp;
110  fPressure = pressure;
111 
112  // Initialize theElementVector allocating one
113  // element corresponding to this material
114  maxNbComponents = fNumberOfComponents = fNumberOfElements = 1;
115  fArrayLength = maxNbComponents;
116  theElementVector = new G4ElementVector();
117 
118  const std::vector<G4String> elmnames =
120  G4String enam, snam;
121  G4int iz = G4lrint(z);
122  if(iz < (G4int)elmnames.size()) {
123  snam = elmnames[iz];
124  enam = snam;
125  } else {
126  enam = "ELM_" + name;
127  snam = name;
128  }
129  theElementVector->push_back(new G4Element(enam, snam, z, a));
130 
131  fMassFractionVector = new G4double[1];
132  fMassFractionVector[0] = 1. ;
133  fMassOfMolecule = a/CLHEP::Avogadro;
134 
135  if (fState == kStateUndefined)
136  {
137  if (fDensity > kGasThreshold) { fState = kStateSolid; }
138  else { fState = kStateGas; }
139  }
140 
141  ComputeDerivedQuantities();
142 }
143 
144 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
145 
146 // Constructor to create a material from a List of constituents
147 // (elements and/or materials) added with AddElement or AddMaterial
148 
150  G4int nComponents,
151  G4State state, G4double temp, G4double pressure)
152  : fName(name)
153 {
154  InitializePointers();
155 
156  if (density < universe_mean_density)
157  {
158  G4cout << "--- Warning from G4Material::G4Material()"
159  << " define a material with density=0 is not allowed. \n"
160  << " The material " << name << " will be constructed with the"
161  << " default minimal density: " << universe_mean_density/(g/cm3)
162  << "g/cm3" << G4endl;
163  density = universe_mean_density;
164  }
165 
166  fDensity = density;
167  fState = state;
168  fTemp = temp;
169  fPressure = pressure;
170 
171  maxNbComponents = nComponents;
172  fArrayLength = maxNbComponents;
173  fNumberOfComponents = fNumberOfElements = 0;
174  theElementVector = new G4ElementVector();
175  theElementVector->reserve(maxNbComponents);
176 
177  if (fState == kStateUndefined)
178  {
179  if (fDensity > kGasThreshold) { fState = kStateSolid; }
180  else { fState = kStateGas; }
181  }
182 }
183 
184 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
185 
186 // Constructor to create a material from base material
187 
189  const G4Material* bmat,
190  G4State state, G4double temp, G4double pressure)
191  : fName(name)
192 {
193  InitializePointers();
194 
195  if (density < universe_mean_density)
196  {
197  G4cout << "--- Warning from G4Material::G4Material()"
198  << " define a material with density=0 is not allowed. \n"
199  << " The material " << name << " will be constructed with the"
200  << " default minimal density: " << universe_mean_density/(g/cm3)
201  << "g/cm3" << G4endl;
202  density = universe_mean_density;
203  }
204 
205  fDensity = density;
206  fState = state;
207  fTemp = temp;
208  fPressure = pressure;
209 
210  fBaseMaterial = bmat;
211  fChemicalFormula = fBaseMaterial->GetChemicalFormula();
212  fMassOfMolecule = fBaseMaterial->GetMassOfMolecule();
213 
214  fNumberOfElements = fBaseMaterial->GetNumberOfElements();
215  maxNbComponents = fNumberOfElements;
216  fNumberOfComponents = fNumberOfElements;
217 
218  CopyPointersOfBaseMaterial();
219 }
220 
221 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
222 
223 // Fake default constructor - sets only member data and allocates memory
224 // for usage restricted to object persistency
225 
227  : fName("")
228 {
229  InitializePointers();
230 }
231 
232 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
233 
235 {
236  // G4cout << "### Destruction of material " << fName << " started" <<G4endl;
237  if(fBaseMaterial == nullptr) {
238  delete theElementVector;
239  delete fSandiaTable;
240  //delete fMaterialPropertiesTable;
241  if(fMassFractionVector) { delete [] fMassFractionVector; }
242  if(fAtomsVector) { delete [] fAtomsVector; }
243  }
244  delete fIonisation;
245  if(VecNbOfAtomsPerVolume) { delete [] VecNbOfAtomsPerVolume; }
246 
247  // Remove this material from theMaterialTable.
248  //
249  theMaterialTable[fIndexInTable] = nullptr;
250 }
251 
252 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
253 
254 void G4Material::InitializePointers()
255 {
256  theElementVector = nullptr;
257  fMassFractionVector = nullptr;
258  fAtomsVector = nullptr;
259  fMaterialPropertiesTable = nullptr;
260 
261  VecNbOfAtomsPerVolume = nullptr;
262  fBaseMaterial = nullptr;
263 
264  fChemicalFormula = "";
265 
266  // initilized data members
267  fDensity = 0.0;
268  fState = kStateUndefined;
269  fTemp = 0.0;
270  fPressure = 0.0;
271  maxNbComponents = 0;
272  fArrayLength = 0;
273  fNumberOfComponents = 0;
274  fNumberOfElements = 0;
275  TotNbOfAtomsPerVolume = 0.0;
276  TotNbOfElectPerVolume = 0.0;
277  fRadlen = 0.0;
278  fNuclInterLen = 0.0;
279  fMassOfMolecule = 0.0;
280 
281  fIonisation = nullptr;
282  fSandiaTable = nullptr;
283 
284  // Store in the static Table of Materials
285  fIndexInTable = theMaterialTable.size();
286  for(size_t i=0; i<fIndexInTable; ++i) {
287  if(theMaterialTable[i]->GetName() == fName) {
288  G4cout << "G4Material WARNING: duplicate name of material "
289  << fName << G4endl;
290  break;
291  }
292  }
293  theMaterialTable.push_back(this);
294 }
295 
296 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
297 
298 void G4Material::ComputeDerivedQuantities()
299 {
300  // Header routine to compute various properties of material.
301  //
302 
303  // Number of atoms per volume (per element), total nb of electrons per volume
304  G4double Zi, Ai;
305  TotNbOfAtomsPerVolume = 0.;
306  if (VecNbOfAtomsPerVolume) { delete [] VecNbOfAtomsPerVolume; }
307  VecNbOfAtomsPerVolume = new G4double[fNumberOfElements];
308  TotNbOfElectPerVolume = 0.;
309  for (G4int i=0; i<fNumberOfElements; ++i) {
310  Zi = (*theElementVector)[i]->GetZ();
311  Ai = (*theElementVector)[i]->GetA();
312  VecNbOfAtomsPerVolume[i] = Avogadro*fDensity*fMassFractionVector[i]/Ai;
313  TotNbOfAtomsPerVolume += VecNbOfAtomsPerVolume[i];
314  TotNbOfElectPerVolume += VecNbOfAtomsPerVolume[i]*Zi;
315  }
316 
317  ComputeRadiationLength();
318  ComputeNuclearInterLength();
319 
320  if (fIonisation) { delete fIonisation; }
321  fIonisation = new G4IonisParamMat(this);
322  if (fSandiaTable) { delete fSandiaTable; }
323  fSandiaTable = new G4SandiaTable(this);
324 }
325 
326 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
327 
328 void G4Material::CopyPointersOfBaseMaterial()
329 {
330  G4double factor = fDensity/fBaseMaterial->GetDensity();
331  TotNbOfAtomsPerVolume = factor*fBaseMaterial->GetTotNbOfAtomsPerVolume();
332  TotNbOfElectPerVolume = factor*fBaseMaterial->GetTotNbOfElectPerVolume();
333 
334  theElementVector =
335  const_cast<G4ElementVector*>(fBaseMaterial->GetElementVector());
336  fMassFractionVector =
337  const_cast<G4double*>(fBaseMaterial->GetFractionVector());
338  fAtomsVector = const_cast<G4int*>(fBaseMaterial->GetAtomsVector());
339 
340  const G4double* v = fBaseMaterial->GetVecNbOfAtomsPerVolume();
341  if (VecNbOfAtomsPerVolume) { delete [] VecNbOfAtomsPerVolume; }
342  VecNbOfAtomsPerVolume = new G4double[fNumberOfElements];
343  for (G4int i=0; i<fNumberOfElements; ++i) {
344  VecNbOfAtomsPerVolume[i] = factor*v[i];
345  }
346  fRadlen = fBaseMaterial->GetRadlen()/factor;
347  fNuclInterLen = fBaseMaterial->GetNuclearInterLength()/factor;
348 
349  if (fIonisation) { delete fIonisation; }
350  fIonisation = new G4IonisParamMat(this);
351 
352  fSandiaTable = fBaseMaterial->GetSandiaTable();
353  fIonisation->SetMeanExcitationEnergy(fBaseMaterial->GetIonisation()->GetMeanExcitationEnergy());
354 
355  fMaterialPropertiesTable = fBaseMaterial->GetMaterialPropertiesTable();
356 }
357 
358 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
359 
360 // AddElement -- composition by atom count
361 
362 void G4Material::AddElement(G4Element* element, G4int nAtoms)
363 {
364  // initialization
365  if ( fNumberOfElements == 0 ) {
366  fAtomsVector = new G4int [fArrayLength];
367  fMassFractionVector = new G4double[fArrayLength];
368  }
369 
370  // filling ...
371  if ( fNumberOfElements < maxNbComponents ) {
372  theElementVector->push_back(element);
373  fAtomsVector[fNumberOfElements] = nAtoms;
374  fNumberOfComponents = ++fNumberOfElements;
375  } else {
376  G4cout << "G4Material::AddElement ERROR for " << fName << " nElement= "
377  << fNumberOfElements << G4endl;
378  G4Exception ("G4Material::AddElement()", "mat031", FatalException,
379  "Attempt to add more than the declared number of elements.");
380  }
381  // filled.
382  if ( fNumberOfElements == maxNbComponents ) {
383  // compute proportion by mass
384  G4int i=0;
385  G4double Amol = 0.;
386  for (i=0; i<fNumberOfElements; ++i) {
387  G4double w = fAtomsVector[i]*(*theElementVector)[i]->GetA();
388  Amol += w;
389  fMassFractionVector[i] = w;
390  }
391  for (i=0; i<fNumberOfElements; ++i) {
392  fMassFractionVector[i] /= Amol;
393  }
394 
395  fMassOfMolecule = Amol/Avogadro;
396  ComputeDerivedQuantities();
397  }
398 }
399 
400 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
401 
402 // AddElement -- composition by fraction of mass
403 
404 void G4Material::AddElement(G4Element* element, G4double fraction)
405 {
406  if(fraction < 0.0 || fraction > 1.0) {
407  G4cout << "G4Material::AddElement ERROR for " << fName << " and "
408  << element->GetName() << " mass fraction= " << fraction
409  << " is wrong " << G4endl;
410  G4Exception ("G4Material::AddElement()", "mat032", FatalException,
411  "Attempt to add element with wrong mass fraction");
412  }
413  // initialization
414  if (fNumberOfComponents == 0) {
415  fMassFractionVector = new G4double[fArrayLength];
416  fAtomsVector = new G4int [fArrayLength];
417  }
418  // filling ...
419  if (fNumberOfComponents < maxNbComponents) {
420  G4int el = 0;
421  // Loop checking, 07-Aug-2015, Vladimir Ivanchenko
422  while ((el<fNumberOfElements)&&(element!=(*theElementVector)[el])) { ++el; }
423  if (el<fNumberOfElements) fMassFractionVector[el] += fraction;
424  else {
425  theElementVector->push_back(element);
426  fMassFractionVector[el] = fraction;
427  ++fNumberOfElements;
428  }
429  ++fNumberOfComponents;
430  } else {
431  G4cout << "G4Material::AddElement ERROR for " << fName << " nElement= "
432  << fNumberOfElements << G4endl;
433  G4Exception ("G4Material::AddElement()", "mat033", FatalException,
434  "Attempt to add more than the declared number of elements.");
435  }
436 
437  // filled.
438  if (fNumberOfComponents == maxNbComponents) {
439 
440  G4int i=0;
441  G4double Zmol(0.), Amol(0.);
442  // check sum of weights -- OK?
443  G4double wtSum(0.0);
444  for (i=0; i<fNumberOfElements; ++i) {
445  wtSum += fMassFractionVector[i];
446  Zmol += fMassFractionVector[i]*(*theElementVector)[i]->GetZ();
447  Amol += fMassFractionVector[i]*(*theElementVector)[i]->GetA();
448  }
449  if (std::fabs(1.-wtSum) > perThousand) {
450  G4cerr << "WARNING !! for " << fName << " sum of fractional masses "
451  << wtSum << " is not 1 - results may be wrong"
452  << G4endl;
453  }
454  for (i=0; i<fNumberOfElements; ++i) {
455  fAtomsVector[i] =
456  G4lrint(fMassFractionVector[i]*Amol/(*theElementVector)[i]->GetA());
457  }
458 
459  ComputeDerivedQuantities();
460  }
461 }
462 
463 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
464 
465 // AddMaterial -- composition by fraction of mass
466 
468 {
469  if(fraction < 0.0 || fraction > 1.0) {
470  G4cout << "G4Material::AddMaterial ERROR for " << fName << " and "
471  << material->GetName() << " mass fraction= " << fraction
472  << " is wrong ";
473  G4Exception ("G4Material::AddMaterial()", "mat034", FatalException,
474  "Attempt to add material with wrong mass fraction");
475  }
476  // initialization
477  if (fNumberOfComponents == 0) {
478  fMassFractionVector = new G4double[fArrayLength];
479  fAtomsVector = new G4int [fArrayLength];
480  }
481 
482  G4int nelm = material->GetNumberOfElements();
483 
484  // arrays should be extended
485  if(nelm > 1) {
486  G4int nold = fArrayLength;
487  fArrayLength += nelm - 1;
488  G4double* v1 = new G4double[fArrayLength];
489  G4int* i1 = new G4int[fArrayLength];
490  for(G4int i=0; i<nold; ++i) {
491  v1[i] = fMassFractionVector[i];
492  i1[i] = fAtomsVector[i];
493  }
494  delete [] fAtomsVector;
495  delete [] fMassFractionVector;
496  fMassFractionVector = v1;
497  fAtomsVector = i1;
498  }
499 
500  // filling ...
501  if (fNumberOfComponents < maxNbComponents) {
502  for (G4int elm=0; elm<nelm; ++elm)
503  {
504  G4Element* element = (*(material->GetElementVector()))[elm];
505  G4int el = 0;
506  // Loop checking, 07-Aug-2015, Vladimir Ivanchenko
507  while ((el<fNumberOfElements)&&(element!=(*theElementVector)[el])) el++;
508  if (el < fNumberOfElements) fMassFractionVector[el] += fraction
509  *(material->GetFractionVector())[elm];
510  else {
511  theElementVector->push_back(element);
512  fMassFractionVector[el] = fraction
513  *(material->GetFractionVector())[elm];
514  ++fNumberOfElements;
515  }
516  }
517  ++fNumberOfComponents;
519  fMatComponents[material] = fraction;
520 
521  } else {
522  G4cout << "G4Material::AddMaterial ERROR for " << fName << " nElement= "
523  << fNumberOfElements << G4endl;
524  G4Exception ("G4Material::AddMaterial()", "mat035", FatalException,
525  "Attempt to add more than the declared number of components.");
526  }
527 
528  // filled.
529  if (fNumberOfComponents == maxNbComponents) {
530  G4int i=0;
531  G4double Zmol(0.), Amol(0.);
532  // check sum of weights -- OK?
533  G4double wtSum(0.0);
534  for (i=0; i<fNumberOfElements; ++i) {
535  wtSum += fMassFractionVector[i];
536  Zmol += fMassFractionVector[i]*(*theElementVector)[i]->GetZ();
537  Amol += fMassFractionVector[i]*(*theElementVector)[i]->GetA();
538  }
539  if (std::fabs(1.-wtSum) > perThousand) {
540  G4cout << "G4Material::AddMaterial WARNING !! for " << fName
541  << " sum of fractional masses "
542  << wtSum << " is not 1 - results may be wrong"
543  << G4endl;
544  }
545  for (i=0; i<fNumberOfElements; ++i) {
546  fAtomsVector[i] =
547  G4lrint(fMassFractionVector[i]*Amol/(*theElementVector)[i]->GetA());
548  }
549 
550  ComputeDerivedQuantities();
551  }
552 }
553 
554 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
555 
556 void G4Material::ComputeRadiationLength()
557 {
558  G4double radinv = 0.0 ;
559  for (G4int i=0;i<fNumberOfElements;++i) {
560  radinv += VecNbOfAtomsPerVolume[i]*((*theElementVector)[i]->GetfRadTsai());
561  }
562  fRadlen = (radinv <= 0.0 ? DBL_MAX : 1./radinv);
563 }
564 
565 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
566 
567 void G4Material::ComputeNuclearInterLength()
568 {
569  static const G4double lambda0 = 35*CLHEP::g/CLHEP::cm2;
570  static const G4double twothird = 2.0/3.0;
571  G4double NILinv = 0.0;
572  for (G4int i=0; i<fNumberOfElements; ++i) {
573  G4int Z = G4lrint( (*theElementVector)[i]->GetZ());
574  G4double A = (*theElementVector)[i]->GetN();
575  if(1 == Z) {
576  NILinv += VecNbOfAtomsPerVolume[i]*A;
577  } else {
578  NILinv += VecNbOfAtomsPerVolume[i]*G4Exp(twothird*G4Log(A));
579  }
580  }
581  NILinv *= amu/lambda0;
582  fNuclInterLen = (NILinv <= 0.0 ? DBL_MAX : 1./NILinv);
583 }
584 
585 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
586 
588 {
589  return &theMaterialTable;
590 }
591 
592 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
593 
595 {
596  return theMaterialTable.size();
597 }
598 
599 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
600 
601 G4Material*
602 G4Material::GetMaterial(const G4String& materialName, G4bool warning)
603 {
604  // search the material by its name
605  for (size_t J=0 ; J<theMaterialTable.size() ; ++J)
606  {
607  if (theMaterialTable[J]->GetName() == materialName)
608  { return theMaterialTable[J]; }
609  }
610 
611  // the material does not exist in the table
612  if (warning) {
613  G4cout << "G4Material::GetMaterial() WARNING: The material: "
614  << materialName
615  << " does not exist in the table. Return NULL pointer."
616  << G4endl;
617  }
618  return nullptr;
619 }
620 
621 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
622 
624 {
625  if (fNumberOfElements > 1) {
626  G4cout << "G4Material ERROR in GetZ. The material: " << fName
627  << " is a mixture.";
628  G4Exception ("G4Material::GetZ()", "mat036", FatalException,
629  "the Atomic number is not well defined." );
630  }
631  return (*theElementVector)[0]->GetZ();
632 }
633 
634 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
635 
637 {
638  if (fNumberOfElements > 1) {
639  G4cout << "G4Material ERROR in GetA. The material: " << fName
640  << " is a mixture.";
641  G4Exception ("G4Material::GetA()", "mat037", FatalException,
642  "the Atomic mass is not well defined." );
643  }
644  return (*theElementVector)[0]->GetA();
645 }
646 
647 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
648 #include "G4ExtendedMaterial.hh"
649 
650 std::ostream& operator<<(std::ostream& flux, const G4Material* material)
651 {
652  std::ios::fmtflags mode = flux.flags();
653  flux.setf(std::ios::fixed,std::ios::floatfield);
654  G4long prec = flux.precision(3);
655 
656  flux
657  << " Material: " << std::setw(8) << material->fName
658  << " " << material->fChemicalFormula << " "
659  << " density: " << std::setw(6) << std::setprecision(3)
660  << G4BestUnit(material->fDensity,"Volumic Mass")
661  << " RadL: " << std::setw(7) << std::setprecision(3)
662  << G4BestUnit(material->fRadlen,"Length")
663  << " Nucl.Int.Length: " << std::setw(7) << std::setprecision(3)
664  << G4BestUnit(material->fNuclInterLen,"Length")
665  << "\n" << std::setw(30)
666  << " Imean: " << std::setw(7) << std::setprecision(3)
668  "Energy");
669 
670  if(material->fState == kStateGas) {
671  flux
672  << " temperature: " << std::setw(6) << std::setprecision(2)
673  << (material->fTemp)/kelvin << " K"
674  << " pressure: " << std::setw(6) << std::setprecision(2)
675  << (material->fPressure)/atmosphere << " atm";
676  }
677  flux << "\n";
678 
679  for (G4int i=0; i<material->fNumberOfElements; i++) {
680  flux
681  << "\n ---> " << (*(material->theElementVector))[i]
682  << "\n ElmMassFraction: "
683  << std::setw(6)<< std::setprecision(2)
684  << (material->fMassFractionVector[i])/perCent << " %"
685  << " ElmAbundance " << std::setw(6)<< std::setprecision(2)
686  << 100*(material->VecNbOfAtomsPerVolume[i])
687  /(material->TotNbOfAtomsPerVolume)
688  << " % \n";
689  }
690  flux.precision(prec);
691  flux.setf(mode,std::ios::floatfield);
692 
693  if(material->IsExtended())
694  { static_cast<const G4ExtendedMaterial*>(material)->Print(flux); }
695 
696  return flux;
697 }
698 
699 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
700 
701 std::ostream& operator<<(std::ostream& flux, const G4Material& material)
702 {
703  flux << &material;
704  return flux;
705 }
706 
707 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
708 
709 std::ostream& operator<<(std::ostream& flux, G4MaterialTable MaterialTable)
710 {
711  //Dump info for all known materials
712  flux << "\n***** Table : Nb of materials = " << MaterialTable.size()
713  << " *****\n" << G4endl;
714 
715  for (size_t i=0; i<MaterialTable.size(); ++i) {
716  flux << MaterialTable[i] << G4endl << G4endl;
717  }
718 
719  return flux;
720 }
721 
722 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
723 
725 { return false; }
726 
727 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
const XML_Char * name
Definition: expat.h:151
G4IonisParamMat * GetIonisation() const
Definition: G4Material.hh:226
G4double GetZ() const
Definition: G4Material.cc:623
G4double GetTotNbOfElectPerVolume() const
Definition: G4Material.hh:212
std::vector< G4Element * > G4ElementVector
G4String fName
Definition: G4AttUtils.hh:55
int kGasThreshold
Definition: hepunit.py:304
static constexpr double cm2
void AddMaterial(G4Material *material, G4double fraction)
Definition: G4Material.cc:467
G4State
Definition: G4Material.hh:114
std::vector< ExP01TrackerHit * > a
Definition: ExP01Classes.hh:33
static G4Material * GetMaterial(const G4String &name, G4bool warning=true)
Definition: G4Material.cc:602
int universe_mean_density
Definition: hepunit.py:307
const G4String & GetChemicalFormula() const
Definition: G4Material.hh:179
void SetMeanExcitationEnergy(G4double value)
static constexpr double perCent
Definition: G4SIunits.hh:332
const G4String & GetName() const
Definition: G4Material.hh:178
const std::vector< G4String > & GetNistElementNames() const
static G4MaterialTable * GetMaterialTable()
Definition: G4Material.cc:587
std::vector< G4Material * > G4MaterialTable
G4double GetDensity() const
Definition: G4Material.hh:180
long G4long
Definition: G4Types.hh:80
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:190
int G4int
Definition: G4Types.hh:78
static constexpr double Avogadro
static G4NistManager * Instance()
float Avogadro
Definition: hepunit.py:253
static constexpr double g
G4SandiaTable * GetSandiaTable() const
Definition: G4Material.hh:229
string material
Definition: eplot.py:19
function g(Y1, Y2, PT2)
Definition: hijing1.383.f:5205
virtual G4bool IsExtended() const
Definition: G4Material.cc:724
const G4double * GetVecNbOfAtomsPerVolume() const
Definition: G4Material.hh:206
static const double prec
Definition: RanecuEngine.cc:58
G4GLOB_DLL std::ostream G4cout
double A(double temperature)
bool G4bool
Definition: G4Types.hh:79
void Print(const std::vector< T > &data)
Definition: DicomRun.hh:109
static size_t GetNumberOfMaterials()
Definition: G4Material.cc:594
static constexpr double kelvin
Definition: G4SIunits.hh:281
G4double GetRadlen() const
Definition: G4Material.hh:220
static constexpr double cm3
Definition: G4SIunits.hh:121
G4Material(const G4String &name, G4double z, G4double a, G4double density, G4State state=kStateUndefined, G4double temp=NTP_Temperature, G4double pressure=CLHEP::STP_Pressure)
Definition: G4Material.cc:90
tuple v
Definition: test.py:18
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *comments)
Definition: G4Exception.cc:41
G4double GetA() const
Definition: G4Material.cc:636
G4double G4Log(G4double x)
Definition: G4Log.hh:230
G4double G4Exp(G4double initial_x)
Exponential Function double precision.
Definition: G4Exp.hh:183
G4double GetTotNbOfAtomsPerVolume() const
Definition: G4Material.hh:209
int G4lrint(double ad)
Definition: templates.hh:163
G4MaterialPropertiesTable * GetMaterialPropertiesTable() const
Definition: G4Material.hh:252
G4double GetMassOfMolecule() const
Definition: G4Material.hh:242
std::ostream & operator<<(std::ostream &, const BasicVector3D< float > &)
const G4int * GetAtomsVector() const
Definition: G4Material.hh:198
tuple z
Definition: test.py:28
G4double GetMeanExcitationEnergy() const
#define G4endl
Definition: G4ios.hh:61
void AddElement(G4Element *element, G4int nAtoms)
Definition: G4Material.cc:362
size_t GetNumberOfElements() const
Definition: G4Material.hh:186
static constexpr double perThousand
Definition: G4SIunits.hh:333
double G4double
Definition: G4Types.hh:76
const G4String & GetName() const
Definition: G4Element.hh:127
const G4double * GetFractionVector() const
Definition: G4Material.hh:194
#define DBL_MAX
Definition: templates.hh:83
virtual ~G4Material()
Definition: G4Material.cc:234
static constexpr double atmosphere
Definition: G4SIunits.hh:237
G4double GetNuclearInterLength() const
Definition: G4Material.hh:223
G4GLOB_DLL std::ostream G4cerr