Geant4_10
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 70847 2013-06-06 11:56:34Z 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 "G4UnitsTable.hh"
78 #include "G4PhysicalConstants.hh"
79 #include "G4SystemOfUnits.hh"
80 
81 G4MaterialTable G4Material::theMaterialTable;
82 
83 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
84 
85 // Constructor to create a material from scratch
86 
89  G4State state, G4double temp, G4double pressure)
90  : fName(name)
91 {
92  InitializePointers();
93 
94  if (density < universe_mean_density)
95  {
96  G4cout << " G4Material WARNING:"
97  << " define a material with density=0 is not allowed. \n"
98  << " The material " << name << " will be constructed with the"
99  << " default minimal density: " << universe_mean_density/(g/cm3)
100  << "g/cm3" << G4endl;
101  density = universe_mean_density;
102  }
103 
104  fDensity = density;
105  fState = state;
106  fTemp = temp;
107  fPressure = pressure;
108 
109  // Initialize theElementVector allocating one
110  // element corresponding to this material
111  maxNbComponents = fNumberOfComponents = fNumberOfElements = 1;
112  fArrayLength = maxNbComponents;
113  fImplicitElement = true;
114  theElementVector = new G4ElementVector();
115  theElementVector->push_back( new G4Element(name, " ", z, a));
116  fMassFractionVector = new G4double[1];
117  fMassFractionVector[0] = 1. ;
118  fMassOfMolecule = a/Avogadro;
119 
120  if (fState == kStateUndefined)
121  {
122  if (fDensity > kGasThreshold) { fState = kStateSolid; }
123  else { fState = kStateGas; }
124  }
125 
126  ComputeDerivedQuantities();
127 }
128 
129 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
130 
131 // Constructor to create a material from a List of constituents
132 // (elements and/or materials) added with AddElement or AddMaterial
133 
136  G4State state, G4double temp, G4double pressure)
137  : fName(name)
138 {
139  InitializePointers();
140 
141  if (density < universe_mean_density)
142  {
143  G4cout << "--- Warning from G4Material::G4Material()"
144  << " define a material with density=0 is not allowed. \n"
145  << " The material " << name << " will be constructed with the"
146  << " default minimal density: " << universe_mean_density/(g/cm3)
147  << "g/cm3" << G4endl;
148  density = universe_mean_density;
149  }
150 
151  fDensity = density;
152  fState = state;
153  fTemp = temp;
154  fPressure = pressure;
155 
156  maxNbComponents = nComponents;
157  fArrayLength = maxNbComponents;
158  fNumberOfComponents = fNumberOfElements = 0;
159  theElementVector = new G4ElementVector();
160  theElementVector->reserve(maxNbComponents);
161 
162  if (fState == kStateUndefined)
163  {
164  if (fDensity > kGasThreshold) { fState = kStateSolid; }
165  else { fState = kStateGas; }
166  }
167 }
168 
169 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
170 
171 // Constructor to create a material from base material
172 
174  const G4Material* bmat,
175  G4State state, G4double temp, G4double pressure)
176  : fName(name)
177 {
178  InitializePointers();
179 
180  if (density < universe_mean_density)
181  {
182  G4cout << "--- Warning from G4Material::G4Material()"
183  << " define a material with density=0 is not allowed. \n"
184  << " The material " << name << " will be constructed with the"
185  << " default minimal density: " << universe_mean_density/(g/cm3)
186  << "g/cm3" << G4endl;
187  density = universe_mean_density;
188  }
189 
190  fDensity = density;
191  fState = state;
192  fTemp = temp;
193  fPressure = pressure;
194 
195  fBaseMaterial = bmat;
196  fChemicalFormula = fBaseMaterial->GetChemicalFormula();
197  fMassOfMolecule = fBaseMaterial->GetMassOfMolecule();
198 
199  fNumberOfElements = fBaseMaterial->GetNumberOfElements();
200  maxNbComponents = fNumberOfElements;
201  fNumberOfComponents = fNumberOfElements;
202 
203  fMaterialPropertiesTable = fBaseMaterial->GetMaterialPropertiesTable();
204 
205  CopyPointersOfBaseMaterial();
206 }
207 
208 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
209 
210 // Fake default constructor - sets only member data and allocates memory
211 // for usage restricted to object persistency
212 
214  : fNumberOfComponents(0), fNumberOfElements(0), theElementVector(0),
215  fImplicitElement(false), fMassFractionVector(0), fAtomsVector(0),
216  fMaterialPropertiesTable(0), fIndexInTable(0),
217  VecNbOfAtomsPerVolume(0)
218 {
219  InitializePointers();
220 }
221 
222 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
223 
225 {
226  // G4cout << "### Destruction of material " << fName << " started" <<G4endl;
227  if(!fBaseMaterial) {
228  if (theElementVector) { delete theElementVector; }
229  if (fMassFractionVector) { delete [] fMassFractionVector; }
230  if (fAtomsVector) { delete [] fAtomsVector; }
231  if (fSandiaTable) { delete fSandiaTable; }
232  }
233  if (fIonisation) { delete fIonisation; }
234  if (VecNbOfAtomsPerVolume) { delete [] VecNbOfAtomsPerVolume; }
235 
236  // Remove this material from theMaterialTable.
237  //
238  theMaterialTable[fIndexInTable] = 0;
239 }
240 
241 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
242 
243 void G4Material::InitializePointers()
244 {
245  theElementVector = 0;
246  fMassFractionVector = 0;
247  fAtomsVector = 0;
248  fMaterialPropertiesTable = 0;
249 
250  VecNbOfAtomsPerVolume = 0;
251  fBaseMaterial = 0;
252 
253  fImplicitElement = false;
254  fChemicalFormula = "";
255 
256  // initilized data members
257  fDensity = 0.0;
258  fState = kStateUndefined;
259  fTemp = 0.0;
260  fPressure = 0.0;
261  maxNbComponents = 0;
262  fArrayLength = 0;
263  TotNbOfAtomsPerVolume = 0;
264  TotNbOfElectPerVolume = 0;
265  fRadlen = 0.0;
266  fNuclInterLen = 0.0;
267  fMassOfMolecule = 0.0;
268 
269  fIonisation = 0;
270  fSandiaTable = 0;
271 
272  // Store in the static Table of Materials
273  fIndexInTable = theMaterialTable.size();
274  for(size_t i=0; i<fIndexInTable; ++i) {
275  if(theMaterialTable[i]->GetName() == fName) {
276  G4cout << "G4Material WARNING: doublicate name of the new material "
277  << fName << G4endl;
278  break;
279  }
280  }
281  theMaterialTable.push_back(this);
282 }
283 
284 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
285 
286 void G4Material::ComputeDerivedQuantities()
287 {
288  // Header routine to compute various properties of material.
289  //
290 
291  // Number of atoms per volume (per element), total nb of electrons per volume
292  G4double Zi, Ai;
293  TotNbOfAtomsPerVolume = 0.;
294  if (VecNbOfAtomsPerVolume) { delete [] VecNbOfAtomsPerVolume; }
295  VecNbOfAtomsPerVolume = new G4double[fNumberOfElements];
296  TotNbOfElectPerVolume = 0.;
297  for (size_t i=0; i<fNumberOfElements; ++i) {
298  Zi = (*theElementVector)[i]->GetZ();
299  Ai = (*theElementVector)[i]->GetA();
300  VecNbOfAtomsPerVolume[i] = Avogadro*fDensity*fMassFractionVector[i]/Ai;
301  TotNbOfAtomsPerVolume += VecNbOfAtomsPerVolume[i];
302  TotNbOfElectPerVolume += VecNbOfAtomsPerVolume[i]*Zi;
303  }
304 
305  ComputeRadiationLength();
306  ComputeNuclearInterLength();
307 
308  if (fIonisation) { delete fIonisation; }
309  fIonisation = new G4IonisParamMat(this);
310  if (fSandiaTable) { delete fSandiaTable; }
311  fSandiaTable = new G4SandiaTable(this);
312 }
313 
314 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
315 
316 void G4Material::CopyPointersOfBaseMaterial()
317 {
318  G4double factor = fDensity/fBaseMaterial->GetDensity();
319  TotNbOfAtomsPerVolume = factor*fBaseMaterial->GetTotNbOfAtomsPerVolume();
320  TotNbOfElectPerVolume = factor*fBaseMaterial->GetTotNbOfElectPerVolume();
321 
322  theElementVector = const_cast<G4ElementVector*>(fBaseMaterial->GetElementVector());
323  fMassFractionVector = const_cast<G4double*>(fBaseMaterial->GetFractionVector());
324  fAtomsVector = const_cast<G4int*>(fBaseMaterial->GetAtomsVector());
325 
326  const G4double* v = fBaseMaterial->GetVecNbOfAtomsPerVolume();
327  if (VecNbOfAtomsPerVolume) { delete [] VecNbOfAtomsPerVolume; }
328  VecNbOfAtomsPerVolume = new G4double[fNumberOfElements];
329  for (size_t i=0; i<fNumberOfElements; ++i) {
330  VecNbOfAtomsPerVolume[i] = factor*v[i];
331  }
332  fRadlen = fBaseMaterial->GetRadlen()/factor;
333  fNuclInterLen = fBaseMaterial->GetNuclearInterLength()/factor;
334  if (fIonisation) { delete fIonisation; }
335  fIonisation = new G4IonisParamMat(this);
336 
337  fSandiaTable = fBaseMaterial->GetSandiaTable();
338  fIonisation->SetMeanExcitationEnergy(fBaseMaterial->GetIonisation()->GetMeanExcitationEnergy());
339 }
340 
341 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
342 
343 // AddElement -- composition by atom count
344 
345 void G4Material::AddElement(G4Element* element, G4int nAtoms)
346 {
347  // initialization
348  if ( fNumberOfElements == 0 ) {
349  fAtomsVector = new G4int [fArrayLength];
350  fMassFractionVector = new G4double[fArrayLength];
351  }
352 
353  // filling ...
354  if ( G4int(fNumberOfElements) < maxNbComponents ) {
355  theElementVector->push_back(element);
356  fAtomsVector[fNumberOfElements] = nAtoms;
357  fNumberOfComponents = ++fNumberOfElements;
358  } else {
359  G4cout << "G4Material::AddElement ERROR for " << fName << " nElement= "
360  << fNumberOfElements << G4endl;
361  G4Exception ("G4Material::AddElement()", "mat031", FatalException,
362  "Attempt to add more than the declared number of elements.");
363  }
364  // filled.
365  if ( G4int(fNumberOfElements) == maxNbComponents ) {
366  // compute proportion by mass
367  size_t i=0;
368  G4double Amol = 0.;
369  for (i=0; i<fNumberOfElements; ++i) {
370  G4double w = fAtomsVector[i]*(*theElementVector)[i]->GetA();
371  Amol += w;
372  fMassFractionVector[i] = w;
373  }
374  for (i=0; i<fNumberOfElements; ++i) {
375  fMassFractionVector[i] /= Amol;
376  }
377 
378  fMassOfMolecule = Amol/Avogadro;
379  ComputeDerivedQuantities();
380  }
381 }
382 
383 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
384 
385 // AddElement -- composition by fraction of mass
386 
387 void G4Material::AddElement(G4Element* element, G4double fraction)
388 {
389  if(fraction < 0.0 || fraction > 1.0) {
390  G4cout << "G4Material::AddElement ERROR for " << fName << " and "
391  << element->GetName() << " mass fraction= " << fraction
392  << " is wrong " << G4endl;
393  G4Exception ("G4Material::AddElement()", "mat032", FatalException,
394  "Attempt to add element with wrong mass fraction");
395  }
396  // initialization
397  if (fNumberOfComponents == 0) {
398  fMassFractionVector = new G4double[fArrayLength];
399  fAtomsVector = new G4int [fArrayLength];
400  }
401  // filling ...
402  if (G4int(fNumberOfComponents) < maxNbComponents) {
403  size_t el = 0;
404  while ((el<fNumberOfElements)&&(element!=(*theElementVector)[el])) { ++el; }
405  if (el<fNumberOfElements) fMassFractionVector[el] += fraction;
406  else {
407  theElementVector->push_back(element);
408  fMassFractionVector[el] = fraction;
409  ++fNumberOfElements;
410  // element->increaseCountUse();
411  }
412  ++fNumberOfComponents;
413  } else {
414  G4cout << "G4Material::AddElement ERROR for " << fName << " nElement= "
415  << fNumberOfElements << G4endl;
416  G4Exception ("G4Material::AddElement()", "mat033", FatalException,
417  "Attempt to add more than the declared number of elements.");
418  }
419 
420  // filled.
421  if (G4int(fNumberOfComponents) == maxNbComponents) {
422 
423  size_t i=0;
424  G4double Zmol(0.), Amol(0.);
425  // check sum of weights -- OK?
426  G4double wtSum(0.0);
427  for (i=0; i<fNumberOfElements; ++i) {
428  wtSum += fMassFractionVector[i];
429  Zmol += fMassFractionVector[i]*(*theElementVector)[i]->GetZ();
430  Amol += fMassFractionVector[i]*(*theElementVector)[i]->GetA();
431  }
432  if (std::fabs(1.-wtSum) > perThousand) {
433  G4cerr << "WARNING !! for " << fName << " sum of fractional masses "
434  << wtSum << " is not 1 - results may be wrong"
435  << G4endl;
436  }
437  for (i=0; i<fNumberOfElements; ++i) {
438  fAtomsVector[i] =
439  G4lrint(fMassFractionVector[i]*Amol/(*theElementVector)[i]->GetA());
440  }
441 
442  ComputeDerivedQuantities();
443  }
444 }
445 
446 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
447 
448 // AddMaterial -- composition by fraction of mass
449 
451 {
452  if(fraction < 0.0 || fraction > 1.0) {
453  G4cout << "G4Material::AddMaterial ERROR for " << fName << " and "
454  << material->GetName() << " mass fraction= " << fraction
455  << " is wrong ";
456  G4Exception ("G4Material::AddMaterial()", "mat034", FatalException,
457  "Attempt to add material with wrong mass fraction");
458  }
459  // initialization
460  if (fNumberOfComponents == 0) {
461  fMassFractionVector = new G4double[fArrayLength];
462  fAtomsVector = new G4int [fArrayLength];
463  }
464 
465  size_t nelm = material->GetNumberOfElements();
466 
467  // arrays should be extended
468  if(nelm > 1) {
469  G4int nold = fArrayLength;
470  fArrayLength += nelm - 1;
471  G4double* v1 = new G4double[fArrayLength];
472  G4int* i1 = new G4int[fArrayLength];
473  for(G4int i=0; i<nold; ++i) {
474  v1[i] = fMassFractionVector[i];
475  i1[i] = fAtomsVector[i];
476  }
477  delete [] fAtomsVector;
478  delete [] fMassFractionVector;
479  fMassFractionVector = v1;
480  fAtomsVector = i1;
481  }
482 
483  // filling ...
484  if (G4int(fNumberOfComponents) < maxNbComponents) {
485  for (size_t elm=0; elm<nelm; ++elm)
486  {
487  G4Element* element = (*(material->GetElementVector()))[elm];
488  size_t el = 0;
489  while ((el<fNumberOfElements)&&(element!=(*theElementVector)[el])) el++;
490  if (el < fNumberOfElements) fMassFractionVector[el] += fraction
491  *(material->GetFractionVector())[elm];
492  else {
493  theElementVector->push_back(element);
494  fMassFractionVector[el] = fraction
495  *(material->GetFractionVector())[elm];
496  ++fNumberOfElements;
497  //element->increaseCountUse();
498  }
499  }
500  ++fNumberOfComponents;
502  fMatComponents[material] = fraction;
503 
504  } else {
505  G4cout << "G4Material::AddMaterial ERROR for " << fName << " nElement= "
506  << fNumberOfElements << G4endl;
507  G4Exception ("G4Material::AddMaterial()", "mat035", FatalException,
508  "Attempt to add more than the declared number of components.");
509  }
510 
511  // filled.
512  if (G4int(fNumberOfComponents) == maxNbComponents) {
513  size_t i=0;
514  G4double Zmol(0.), Amol(0.);
515  // check sum of weights -- OK?
516  G4double wtSum(0.0);
517  for (i=0; i<fNumberOfElements; ++i) {
518  wtSum += fMassFractionVector[i];
519  Zmol += fMassFractionVector[i]*(*theElementVector)[i]->GetZ();
520  Amol += fMassFractionVector[i]*(*theElementVector)[i]->GetA();
521  }
522  if (std::fabs(1.-wtSum) > perThousand) {
523  G4cout << "G4Material::AddMaterial WARNING !! for " << fName
524  << " sum of fractional masses "
525  << wtSum << " is not 1 - results may be wrong"
526  << G4endl;
527  }
528  for (i=0;i<fNumberOfElements;i++) {
529  fAtomsVector[i] =
530  G4lrint(fMassFractionVector[i]*Amol/(*theElementVector)[i]->GetA());
531  }
532 
533  ComputeDerivedQuantities();
534  }
535 }
536 
537 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
538 
539 void G4Material::ComputeRadiationLength()
540 {
541  G4double radinv = 0.0 ;
542  for (size_t i=0;i<fNumberOfElements;++i) {
543  radinv += VecNbOfAtomsPerVolume[i]*((*theElementVector)[i]->GetfRadTsai());
544  }
545  fRadlen = (radinv <= 0.0 ? DBL_MAX : 1./radinv);
546 }
547 
548 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
549 
550 void G4Material::ComputeNuclearInterLength()
551 {
552  static const G4double lambda0 = 35*g/cm2;
553  G4double NILinv = 0.0;
554  for (size_t i=0; i<fNumberOfElements; ++i) {
555  NILinv +=
556  VecNbOfAtomsPerVolume[i]*std::pow((*theElementVector)[i]->GetN(),0.6666666667);
557  }
558  NILinv *= amu/lambda0;
559  fNuclInterLen = (NILinv <= 0.0 ? DBL_MAX : 1./NILinv);
560 }
561 
562 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
563 
565 {
566  return &theMaterialTable;
567 }
568 
569 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
570 
572 {
573  return theMaterialTable.size();
574 }
575 
576 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
577 
578 G4Material* G4Material::GetMaterial(const G4String& materialName, G4bool warning)
579 {
580  // search the material by its name
581  for (size_t J=0 ; J<theMaterialTable.size() ; ++J)
582  {
583  if (theMaterialTable[J]->GetName() == materialName)
584  { return theMaterialTable[J]; }
585  }
586 
587  // the material does not exist in the table
588  if (warning) {
589  G4cout << "G4Material::GetMaterial() WARNING: The material: "
590  << materialName << " does not exist in the table. Return NULL pointer."
591  << G4endl;
592  }
593  return 0;
594 }
595 
596 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
597 
599 {
600  InitializePointers();
601  *this = right;
602 }
603 
604 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
605 
607 {
608  if (fNumberOfElements > 1) {
609  G4cout << "G4Material ERROR in GetZ. The material: " << fName
610  << " is a mixture.";
611  G4Exception ("G4Material::GetZ()", "mat036", FatalException,
612  "the Atomic number is not well defined." );
613  }
614  return (*theElementVector)[0]->GetZ();
615 }
616 
617 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
618 
620 {
621  if (fNumberOfElements > 1) {
622  G4cout << "G4Material ERROR in GetA. The material: " << fName
623  << " is a mixture.";
624  G4Exception ("G4Material::GetA()", "mat037", FatalException,
625  "the Atomic mass is not well defined." );
626  }
627  return (*theElementVector)[0]->GetA();
628 }
629 
630 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
631 
632 const G4Material& G4Material::operator=(const G4Material& right)
633 {
634  if (this != &right)
635  {
636  fName = right.fName;
637  fChemicalFormula = right.fChemicalFormula;
638  fDensity = right.fDensity;
639  fState = right.fState;
640  fTemp = right.fTemp;
641  fPressure = right.fPressure;
642 
643  if(!fBaseMaterial) {
644  if (theElementVector) { delete theElementVector; }
645  if (fMassFractionVector) { delete [] fMassFractionVector; }
646  if (fAtomsVector) { delete [] fAtomsVector; }
647  if (fIonisation) { delete fIonisation; }
648  if (fSandiaTable) { delete fSandiaTable; }
649  }
650 
651  if (VecNbOfAtomsPerVolume) { delete [] VecNbOfAtomsPerVolume; }
652 
653  maxNbComponents = right.maxNbComponents;
654  fNumberOfComponents = right.fNumberOfComponents;
655  fNumberOfElements = right.fNumberOfElements;
656  fImplicitElement = right.fImplicitElement;
657 
658  fMaterialPropertiesTable = right.fMaterialPropertiesTable;
659  fBaseMaterial = right.fBaseMaterial;
660  fMassOfMolecule= right.fMassOfMolecule;
661  fMatComponents= right.fMatComponents;
662 
663  if(fBaseMaterial) {
664  CopyPointersOfBaseMaterial();
665 
666  } else {
667  theElementVector = new G4ElementVector(fNumberOfElements,0);
668  fMassFractionVector = new G4double[fNumberOfElements];
669  fAtomsVector = new G4int[fNumberOfElements];
670  for (size_t i=0; i<fNumberOfElements; ++i) {
671  (*theElementVector)[i] = (*right.theElementVector)[i];
672  fMassFractionVector[i] = right.fMassFractionVector[i];
673  fAtomsVector[i] = right.fAtomsVector[i];
674  }
675  ComputeDerivedQuantities();
676  }
677  }
678  return *this;
679 }
680 
681 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
682 
684 {
685  return (this == (G4Material *) &right);
686 }
687 
688 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
689 
691 {
692  return (this != (G4Material *) &right);
693 }
694 
695 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
696 
697 std::ostream& operator<<(std::ostream& flux, G4Material* material)
698 {
699  std::ios::fmtflags mode = flux.flags();
700  flux.setf(std::ios::fixed,std::ios::floatfield);
701  G4long prec = flux.precision(3);
702 
703  flux
704  << " Material: " << std::setw(8) << material->fName
705  << " " << material->fChemicalFormula << " "
706  << " density: " << std::setw(6) << std::setprecision(3)
707  << G4BestUnit(material->fDensity,"Volumic Mass")
708  << " RadL: " << std::setw(7) << std::setprecision(3)
709  << G4BestUnit(material->fRadlen,"Length")
710  << " Nucl.Int.Length: " << std::setw(7) << std::setprecision(3)
711  << G4BestUnit(material->fNuclInterLen,"Length") <<"\n" << std::setw(30)
712  << " Imean: " << std::setw(7) << std::setprecision(3)
713  << G4BestUnit(material->GetIonisation()->GetMeanExcitationEnergy(),"Energy");
714 
715  if(material->fState == kStateGas) {
716  flux
717  << " temperature: " << std::setw(6) << std::setprecision(2)
718  << (material->fTemp)/kelvin << " K"
719  << " pressure: " << std::setw(6) << std::setprecision(2)
720  << (material->fPressure)/atmosphere << " atm";
721  }
722  flux << "\n";
723 
724  for (size_t i=0; i<material->fNumberOfElements; i++) {
725  flux
726  << "\n ---> " << (*(material->theElementVector))[i]
727  << "\n ElmMassFraction: "
728  << std::setw(6)<< std::setprecision(2)
729  << (material->fMassFractionVector[i])/perCent << " %"
730  << " ElmAbundance " << std::setw(6)<< std::setprecision(2)
731  << 100*(material->VecNbOfAtomsPerVolume[i])/(material->TotNbOfAtomsPerVolume)
732  << " % \n";
733  }
734  flux.precision(prec);
735  flux.setf(mode,std::ios::floatfield);
736 
737  return flux;
738 }
739 
740 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
741 
742 std::ostream& operator<<(std::ostream& flux, G4Material& material)
743 {
744  flux << &material;
745  return flux;
746 }
747 
748 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
749 
750 std::ostream& operator<<(std::ostream& flux, G4MaterialTable MaterialTable)
751 {
752  //Dump info for all known materials
753  flux << "\n***** Table : Nb of materials = " << MaterialTable.size()
754  << " *****\n" << G4endl;
755 
756  for (size_t i=0; i<MaterialTable.size(); ++i) {
757  flux << MaterialTable[i] << G4endl << G4endl;
758  }
759 
760  return flux;
761 }
762 
763 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
G4IonisParamMat * GetIonisation() const
Definition: G4Material.hh:224
G4double GetZ() const
Definition: G4Material.cc:606
G4int operator==(const G4Material &) const
Definition: G4Material.cc:683
G4double GetTotNbOfElectPerVolume() const
Definition: G4Material.hh:210
std::vector< G4Element * > G4ElementVector
tuple a
Definition: test.py:11
G4String fName
Definition: G4AttUtils.hh:55
int kGasThreshold
Definition: hepunit.py:304
void AddMaterial(G4Material *material, G4double fraction)
Definition: G4Material.cc:450
G4State
Definition: G4Material.hh:114
static G4Material * GetMaterial(const G4String &name, G4bool warning=true)
Definition: G4Material.cc:578
int universe_mean_density
Definition: hepunit.py:307
const G4String & GetChemicalFormula() const
Definition: G4Material.hh:177
void SetMeanExcitationEnergy(G4double value)
const G4String & GetName() const
Definition: G4Material.hh:176
static G4MaterialTable * GetMaterialTable()
Definition: G4Material.cc:564
std::vector< G4Material * > G4MaterialTable
G4double GetDensity() const
Definition: G4Material.hh:178
int atmosphere
Definition: hepunit.py:151
long G4long
Definition: G4Types.hh:80
const XML_Char * name
Definition: expat.h:151
float perThousand
Definition: hepunit.py:240
#define G4BestUnit(a, b)
#define G4_USE_G4BESTUNIT_FOR_VERBOSE 1
G4int nComponents
Definition: TRTMaterials.hh:41
const G4ElementVector * GetElementVector() const
Definition: G4Material.hh:188
int G4int
Definition: G4Types.hh:78
float Avogadro
Definition: hepunit.py:253
G4Material(const G4String &name, G4double z, G4double a, G4double density, G4State state=kStateUndefined, G4double temp=CLHEP::STP_Temperature, G4double pressure=CLHEP::STP_Pressure)
Definition: G4Material.cc:87
G4SandiaTable * GetSandiaTable() const
Definition: G4Material.hh:227
string material
Definition: eplot.py:19
G4double density
Definition: TRTMaterials.hh:39
function g(Y1, Y2, PT2)
Definition: hijing1.383.f:5205
const G4double * GetVecNbOfAtomsPerVolume() const
Definition: G4Material.hh:204
G4GLOB_DLL std::ostream G4cout
bool G4bool
Definition: G4Types.hh:79
static size_t GetNumberOfMaterials()
Definition: G4Material.cc:571
G4double GetRadlen() const
Definition: G4Material.hh:218
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:619
G4double GetTotNbOfAtomsPerVolume() const
Definition: G4Material.hh:207
int G4lrint(double ad)
Definition: templates.hh:163
G4int operator!=(const G4Material &) const
Definition: G4Material.cc:690
G4MaterialPropertiesTable * GetMaterialPropertiesTable() const
Definition: G4Material.hh:250
G4double GetMassOfMolecule() const
Definition: G4Material.hh:240
std::ostream & operator<<(std::ostream &, const BasicVector3D< float > &)
const G4int * GetAtomsVector() const
Definition: G4Material.hh:196
tuple z
Definition: test.py:28
float perCent
Definition: hepunit.py:239
G4double GetMeanExcitationEnergy() const
#define G4endl
Definition: G4ios.hh:61
void AddElement(G4Element *element, G4int nAtoms)
Definition: G4Material.cc:345
size_t GetNumberOfElements() const
Definition: G4Material.hh:184
double G4double
Definition: G4Types.hh:76
const G4String & GetName() const
Definition: G4Element.hh:127
const G4double * GetFractionVector() const
Definition: G4Material.hh:192
#define DBL_MAX
Definition: templates.hh:83
virtual ~G4Material()
Definition: G4Material.cc:224
G4double GetNuclearInterLength() const
Definition: G4Material.hh:221
G4GLOB_DLL std::ostream G4cerr