| Geant4 Cross Reference |
1 // 1 //
2 // ******************************************* 2 // ********************************************************************
3 // * License and Disclaimer 3 // * License and Disclaimer *
4 // * 4 // * *
5 // * The Geant4 software is copyright of th 5 // * The Geant4 software is copyright of the Copyright Holders of *
6 // * the Geant4 Collaboration. It is provided 6 // * the Geant4 Collaboration. It is provided under the terms and *
7 // * conditions of the Geant4 Software License 7 // * conditions of the Geant4 Software License, included in the file *
8 // * LICENSE and available at http://cern.ch/ 8 // * LICENSE and available at http://cern.ch/geant4/license . These *
9 // * include a list of copyright holders. 9 // * include a list of copyright holders. *
10 // * 10 // * *
11 // * Neither the authors of this software syst 11 // * Neither the authors of this software system, nor their employing *
12 // * institutes,nor the agencies providing fin 12 // * institutes,nor the agencies providing financial support for this *
13 // * work make any representation or warran 13 // * work make any representation or warranty, express or implied, *
14 // * regarding this software system or assum 14 // * regarding this software system or assume any liability for its *
15 // * use. Please see the license in the file 15 // * use. Please see the license in the file LICENSE and URL above *
16 // * for the full disclaimer and the limitatio 16 // * for the full disclaimer and the limitation of liability. *
17 // * 17 // * *
18 // * This code implementation is the result 18 // * This code implementation is the result of the scientific and *
19 // * technical work of the GEANT4 collaboratio 19 // * technical work of the GEANT4 collaboration. *
20 // * By using, copying, modifying or distri 20 // * By using, copying, modifying or distributing the software (or *
21 // * any work based on the software) you ag 21 // * any work based on the software) you agree to acknowledge its *
22 // * use in resulting scientific publicati 22 // * use in resulting scientific publications, and indicate your *
23 // * acceptance of all terms of the Geant4 Sof 23 // * acceptance of all terms of the Geant4 Software license. *
24 // ******************************************* 24 // ********************************************************************
25 << 25 //
>> 26 // $Id: G4Material.cc 67044 2013-01-30 08:50:06Z gcosmo $
>> 27 //
>> 28 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
>> 29 //
26 // 26-06-96, Code uses operators (+=, *=, ++, 30 // 26-06-96, Code uses operators (+=, *=, ++, -> etc.) correctly, P. Urban
27 // 10-07-96, new data members added by L.Urban 31 // 10-07-96, new data members added by L.Urban
28 // 12-12-96, new data memberfFreeElecDensitys << 32 // 12-12-96, new data members added by L.Urban
29 // 20-01-97, aesthetic rearrangement. RadLengt 33 // 20-01-97, aesthetic rearrangement. RadLength calculation modified.
30 // Data members Zeff and Aeff REMOVE 34 // Data members Zeff and Aeff REMOVED (i.e. passed to the Elements).
31 // (local definition of Zeff in Dens 35 // (local definition of Zeff in DensityEffect and FluctModel...)
32 // Vacuum defined as a G4State. Mixt << 36 // Vacuum defined as a G4State. Mixture flag removed, M.Maire.
33 // 29-01-97, State=Vacuum automatically set de 37 // 29-01-97, State=Vacuum automatically set density=0 in the contructors.
34 // Subsequent protections have been << 38 // Subsequent protections have been put in the calculation of
35 // MeanExcEnergy, ShellCorrectionVec 39 // MeanExcEnergy, ShellCorrectionVector, DensityEffect, M.Maire.
36 // 11-02-97, ComputeDensityEffect() rearranged 40 // 11-02-97, ComputeDensityEffect() rearranged, M.Maire.
37 // 20-03-97, corrected initialization of point 41 // 20-03-97, corrected initialization of pointers, M.Maire.
38 // 28-05-98, the kState=kVacuum has been remov 42 // 28-05-98, the kState=kVacuum has been removed.
39 // automatic check for a minimal den << 43 // automatic check for a minimal density, M.Maire
40 // 12-06-98, new method AddMaterial() allowing << 44 // 12-06-98, new method AddMaterial() allowing mixture of materials, M.Maire
41 // 09-07-98, ionisation parameters removed fro 45 // 09-07-98, ionisation parameters removed from the class, M.Maire
42 // 05-10-98, change names: NumDensity -> NbOfA 46 // 05-10-98, change names: NumDensity -> NbOfAtomsPerVolume
43 // 18-11-98, new interface to SandiaTable 47 // 18-11-98, new interface to SandiaTable
44 // 19-01-99 enlarge tolerance on test of cohe 48 // 19-01-99 enlarge tolerance on test of coherence of gas conditions
45 // 19-07-99, Constructors with chemicalFormula 49 // 19-07-99, Constructors with chemicalFormula added by V.Ivanchenko
46 // 16-01-01, Nuclear interaction length, M.Mai 50 // 16-01-01, Nuclear interaction length, M.Maire
47 // 12-03-01, G4bool fImplicitElement; 51 // 12-03-01, G4bool fImplicitElement;
48 // copy constructor and assignement 52 // copy constructor and assignement operator revised (mma)
49 // 03-05-01, flux.precision(prec) at begin/end 53 // 03-05-01, flux.precision(prec) at begin/end of operator<<
50 // 17-07-01, migration to STL. M. Verderi. 54 // 17-07-01, migration to STL. M. Verderi.
51 // 14-09-01, Suppression of the data member fI 55 // 14-09-01, Suppression of the data member fIndexInTable
52 // 26-02-02, fIndexInTable renewed 56 // 26-02-02, fIndexInTable renewed
53 // 16-04-02, G4Exception put in constructor wi 57 // 16-04-02, G4Exception put in constructor with chemical formula
54 // 06-05-02, remove the check of the ideal gas 58 // 06-05-02, remove the check of the ideal gas state equation
55 // 06-08-02, remove constructors with chemical 59 // 06-08-02, remove constructors with chemical formula (mma)
56 // 22-01-04, proper STL handling of theElement 60 // 22-01-04, proper STL handling of theElementVector (Hisaya)
57 // 30-03-05, warning in GetMaterial(materialNa << 61 // 30-03-05, warning in GetMaterial(materialName)
58 // 09-03-06, minor change of printout (V.Ivanc << 62 // 09-03-06, minor change of printout (V.Ivanchenko)
59 // 10-01-07, compute fAtomVector in the case o << 63 // 10-01-07, compute fAtomVector in the case of mass fraction (V.Ivanchenko)
60 // 27-07-07, improve destructor (V.Ivanchenko) << 64 // 27-07-07, improve destructor (V.Ivanchenko)
61 // 18-10-07, moved definition of mat index to << 65 // 18-10-07, move definition of material index to InitialisePointers (V.Ivanchenko)
62 // 13-08-08, do not use fixed size arrays (V.I 66 // 13-08-08, do not use fixed size arrays (V.Ivanchenko)
63 // 26-10-11, new scheme for G4Exception (mma) 67 // 26-10-11, new scheme for G4Exception (mma)
64 // 13-04-12, map<G4Material*,G4double> fMatCom 68 // 13-04-12, map<G4Material*,G4double> fMatComponents, filled in AddMaterial()
65 // 21-04-12, fMassOfMolecule, computed for Ato 69 // 21-04-12, fMassOfMolecule, computed for AtomsCount (mma)
>> 70 //
>> 71 //
>> 72 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
66 73
67 #include "G4Material.hh" << 74 #include <iomanip>
68 75
69 #include "G4ApplicationState.hh" << 76 #include "G4Material.hh"
70 #include "G4AtomicShells.hh" << 77 #include "G4UnitsTable.hh"
71 #include "G4ExtendedMaterial.hh" <<
72 #include "G4Pow.hh" 78 #include "G4Pow.hh"
73 #include "G4NistManager.hh" <<
74 #include "G4PhysicalConstants.hh" 79 #include "G4PhysicalConstants.hh"
75 #include "G4StateManager.hh" <<
76 #include "G4SystemOfUnits.hh" 80 #include "G4SystemOfUnits.hh"
77 #include "G4UnitsTable.hh" <<
78 81
79 #include <iomanip> << 82 G4MaterialTable G4Material::theMaterialTable;
80 83
81 //....oooOO0OOooo........oooOO0OOooo........oo 84 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
82 85
83 // Constructor to create a material from scrat 86 // Constructor to create a material from scratch
84 87
85 G4Material::G4Material(const G4String& name, G << 88 G4Material::G4Material(const G4String& name, G4double z,
86 G4State state, G4double temp, G4double press << 89 G4double a, G4double density,
87 : fName(name) << 90 G4State state, G4double temp, G4double pressure)
>> 91 : fName(name)
88 { 92 {
89 InitializePointers(); 93 InitializePointers();
90 << 94
91 if (density < universe_mean_density) { << 95 if (density < universe_mean_density)
92 G4cout << " G4Material WARNING:" << 96 {
93 << " define a material with density << 97 G4cout << "--- Warning from G4Material::G4Material()"
94 << " The material " << name << " wi << 98 << " define a material with density=0 is not allowed. \n"
95 << " default minimal density: " << << 99 << " The material " << name << " will be constructed with the"
96 << G4endl; << 100 << " default minimal density: " << universe_mean_density/(g/cm3)
97 density = universe_mean_density; << 101 << "g/cm3" << G4endl;
98 } << 102 density = universe_mean_density;
99 << 103 }
100 fDensity = density; << 104
101 fState = state; << 105 fDensity = density;
102 fTemp = temp; << 106 fState = state;
>> 107 fTemp = temp;
103 fPressure = pressure; 108 fPressure = pressure;
104 109
105 // Initialize theElementVector allocating on 110 // Initialize theElementVector allocating one
106 // element corresponding to this material 111 // element corresponding to this material
107 fNbComponents = fNumberOfElements = 1; << 112 maxNbComponents = fNumberOfComponents = fNumberOfElements = 1;
108 theElementVector = new G4ElementVector(); << 113 fArrayLength = maxNbComponents;
109 << 114 fImplicitElement = true;
110 // take element from DB << 115 theElementVector = new G4ElementVector();
111 G4NistManager* nist = G4NistManager::Instanc << 116 theElementVector->push_back( new G4Element(name, " ", z, a));
112 G4int iz = G4lrint(z); << 117 fMassFractionVector = new G4double[1];
113 auto elm = nist->FindOrBuildElement(iz); << 118 fMassFractionVector[0] = 1. ;
114 if (elm == nullptr) { << 119 fMassOfMolecule = a/Avogadro;
115 elm = new G4Element("ELM_" + name, name, z << 120
116 } << 121 (*theElementVector)[0] -> increaseCountUse();
117 theElementVector->push_back(elm); << 122
118 << 123 if (fState == kStateUndefined)
119 fMassFractionVector = new G4double[1]; << 124 {
120 fMassFractionVector[0] = 1.; << 125 if (fDensity > kGasThreshold) { fState = kStateSolid; }
121 fMassOfMolecule = a / CLHEP::Avogadro; << 126 else { fState = kStateGas; }
122 <<
123 if (fState == kStateUndefined) { <<
124 if (fDensity > kGasThreshold) { <<
125 fState = kStateSolid; <<
126 } <<
127 else { <<
128 fState = kStateGas; <<
129 } 127 }
130 } <<
131 128
132 ComputeDerivedQuantities(); 129 ComputeDerivedQuantities();
133 } 130 }
134 131
135 //....oooOO0OOooo........oooOO0OOooo........oo 132 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
136 133
137 // Constructor to create a material from a Lis 134 // Constructor to create a material from a List of constituents
138 // (elements and/or materials) added with Add 135 // (elements and/or materials) added with AddElement or AddMaterial
139 136
140 G4Material::G4Material(const G4String& name, G << 137 G4Material::G4Material(const G4String& name, G4double density,
141 G4double temp, G4double pressure) << 138 G4int nComponents,
142 : fName(name) << 139 G4State state, G4double temp, G4double pressure)
>> 140 : fName(name)
143 { 141 {
144 InitializePointers(); 142 InitializePointers();
145 << 143
146 if (density < universe_mean_density) { << 144 if (density < universe_mean_density)
147 G4cout << "--- Warning from G4Material::G4 << 145 {
148 << " define a material with density << 146 G4cout << "--- Warning from G4Material::G4Material()"
149 << " The material " << name << " wi << 147 << " define a material with density=0 is not allowed. \n"
150 << " default minimal density: " << << 148 << " The material " << name << " will be constructed with the"
151 << G4endl; << 149 << " default minimal density: " << universe_mean_density/(g/cm3)
152 density = universe_mean_density; << 150 << "g/cm3" << G4endl;
153 } << 151 density = universe_mean_density;
154 << 152 }
155 fDensity = density; << 153
156 fState = state; << 154 fDensity = density;
157 fTemp = temp; << 155 fState = state;
>> 156 fTemp = temp;
158 fPressure = pressure; 157 fPressure = pressure;
159 << 158
160 fNbComponents = nComponents; << 159 maxNbComponents = nComponents;
161 fMassFraction = true; << 160 fArrayLength = maxNbComponents;
162 << 161 fNumberOfComponents = fNumberOfElements = 0;
163 if (fState == kStateUndefined) { << 162 theElementVector = new G4ElementVector();
164 if (fDensity > kGasThreshold) { << 163 theElementVector->reserve(maxNbComponents);
165 fState = kStateSolid; << 164
166 } << 165 if (fState == kStateUndefined)
167 else { << 166 {
168 fState = kStateGas; << 167 if (fDensity > kGasThreshold) { fState = kStateSolid; }
>> 168 else { fState = kStateGas; }
169 } 169 }
170 } <<
171 } 170 }
172 171
173 //....oooOO0OOooo........oooOO0OOooo........oo 172 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
174 173
175 // Constructor to create a material from base 174 // Constructor to create a material from base material
176 175
177 G4Material::G4Material(const G4String& name, G << 176 G4Material::G4Material(const G4String& name, G4double density,
178 G4State state, G4double temp, G4double press << 177 const G4Material* bmat,
179 : fName(name) << 178 G4State state, G4double temp, G4double pressure)
>> 179 : fName(name)
180 { 180 {
181 InitializePointers(); 181 InitializePointers();
>> 182
>> 183 if (density < universe_mean_density)
>> 184 {
>> 185 G4cout << "--- Warning from G4Material::G4Material()"
>> 186 << " define a material with density=0 is not allowed. \n"
>> 187 << " The material " << name << " will be constructed with the"
>> 188 << " default minimal density: " << universe_mean_density/(g/cm3)
>> 189 << "g/cm3" << G4endl;
>> 190 density = universe_mean_density;
>> 191 }
182 192
183 if (density < universe_mean_density) { << 193 fDensity = density;
184 G4cout << "--- Warning from G4Material::G4 << 194 fState = state;
185 << " define a material with density << 195 fTemp = temp;
186 << " The material " << name << " wi <<
187 << " default minimal density: " << <<
188 << G4endl; <<
189 density = universe_mean_density; <<
190 } <<
191 <<
192 fDensity = density; <<
193 fState = state; <<
194 fTemp = temp; <<
195 fPressure = pressure; 196 fPressure = pressure;
196 197
197 fBaseMaterial = bmat; 198 fBaseMaterial = bmat;
198 auto ptr = bmat; <<
199 if (nullptr != ptr) { <<
200 while (true) { <<
201 ptr = ptr->GetBaseMaterial(); <<
202 if (nullptr == ptr) { <<
203 break; <<
204 } <<
205 fBaseMaterial = ptr; <<
206 } <<
207 } <<
208 <<
209 fChemicalFormula = fBaseMaterial->GetChemica 199 fChemicalFormula = fBaseMaterial->GetChemicalFormula();
210 fMassOfMolecule = fBaseMaterial->GetMassOfMo << 200 fMassOfMolecule = fBaseMaterial->GetMassOfMolecule();
211 201
212 fNumberOfElements = (G4int)fBaseMaterial->Ge << 202 fNumberOfElements = fBaseMaterial->GetNumberOfElements();
213 fNbComponents = fNumberOfElements; << 203 maxNbComponents = fNumberOfElements;
>> 204 fNumberOfComponents = fNumberOfElements;
>> 205
>> 206 fMaterialPropertiesTable = fBaseMaterial->GetMaterialPropertiesTable();
214 207
215 CopyPointersOfBaseMaterial(); 208 CopyPointersOfBaseMaterial();
216 } 209 }
217 210
218 //....oooOO0OOooo........oooOO0OOooo........oo 211 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
219 212
>> 213 // Fake default constructor - sets only member data and allocates memory
>> 214 // for usage restricted to object persistency
>> 215
>> 216 G4Material::G4Material(__void__&)
>> 217 : fNumberOfComponents(0), fNumberOfElements(0), theElementVector(0),
>> 218 fImplicitElement(false), fMassFractionVector(0), fAtomsVector(0),
>> 219 fMaterialPropertiesTable(0), fIndexInTable(0),
>> 220 VecNbOfAtomsPerVolume(0), fIonisation(0), fSandiaTable(0)
>> 221 {
>> 222 InitializePointers();
>> 223 }
>> 224
>> 225 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
>> 226
220 G4Material::~G4Material() 227 G4Material::~G4Material()
221 { 228 {
222 if (fBaseMaterial == nullptr) { << 229 // G4cout << "### Destruction of material " << fName << " started" <<G4endl;
223 delete theElementVector; << 230 if(!fBaseMaterial) {
224 delete fSandiaTable; << 231 if (theElementVector) { delete theElementVector; }
225 delete[] fMassFractionVector; << 232 if (fMassFractionVector) { delete [] fMassFractionVector; }
226 delete[] fAtomsVector; << 233 if (fAtomsVector) { delete [] fAtomsVector; }
>> 234 if (fSandiaTable) { delete fSandiaTable; }
227 } 235 }
228 delete fIonisation; << 236 if (fIonisation) { delete fIonisation; }
229 delete[] fVecNbOfAtomsPerVolume; << 237 if (VecNbOfAtomsPerVolume) { delete [] VecNbOfAtomsPerVolume; }
230 238
231 // Remove this material from the MaterialTab << 239 // Remove this material from theMaterialTable.
232 // 240 //
233 (*GetMaterialTable())[fIndexInTable] = nullp << 241 theMaterialTable[fIndexInTable] = 0;
234 } 242 }
235 243
236 //....oooOO0OOooo........oooOO0OOooo........oo 244 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
237 245
238 void G4Material::InitializePointers() 246 void G4Material::InitializePointers()
239 { 247 {
240 fBaseMaterial = nullptr; << 248 theElementVector = 0;
241 fMaterialPropertiesTable = nullptr; << 249 fMassFractionVector = 0;
242 theElementVector = nullptr; << 250 fAtomsVector = 0;
243 fAtomsVector = nullptr; << 251 fMaterialPropertiesTable = 0;
244 fMassFractionVector = nullptr; << 252
245 fVecNbOfAtomsPerVolume = nullptr; << 253 VecNbOfAtomsPerVolume = 0;
246 << 254 fIonisation = 0;
247 fIonisation = nullptr; << 255 fSandiaTable = 0;
248 fSandiaTable = nullptr; << 256
249 << 257 fBaseMaterial = 0;
250 fDensity = fFreeElecDensity = fTemp = fPress << 258
251 fTotNbOfAtomsPerVolume = 0.0; << 259 fImplicitElement = false;
252 fTotNbOfElectPerVolume = 0.0; << 260 fChemicalFormula = "";
253 fRadlen = fNuclInterLen = fMassOfMolecule = << 261
254 << 262 // initilized data members
255 fState = kStateUndefined; << 263 fDensity = 0.0;
256 fNumberOfElements = fNbComponents = fIdxComp << 264 fState = kStateUndefined;
257 fMassFraction = true; << 265 fTemp = 0.0;
258 fChemicalFormula = ""; << 266 fPressure = 0.0;
>> 267 maxNbComponents = 0;
>> 268 fArrayLength = 0;
>> 269 TotNbOfAtomsPerVolume = 0;
>> 270 TotNbOfElectPerVolume = 0;
>> 271 fRadlen = 0.0;
>> 272 fNuclInterLen = 0.0;
>> 273 fMassOfMolecule = 0.0;
259 274
260 // Store in the static Table of Materials 275 // Store in the static Table of Materials
261 fIndexInTable = GetMaterialTable()->size(); << 276 theMaterialTable.push_back(this);
262 for (std::size_t i = 0; i < fIndexInTable; + << 277 fIndexInTable = theMaterialTable.size() - 1;
263 if ((*GetMaterialTable())[i]->GetName() == <<
264 G4cout << "G4Material WARNING: duplicate <<
265 break; <<
266 } <<
267 } <<
268 GetMaterialTable()->push_back(this); <<
269 } 278 }
270 279
271 //....oooOO0OOooo........oooOO0OOooo........oo 280 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
272 281
273 void G4Material::ComputeDerivedQuantities() 282 void G4Material::ComputeDerivedQuantities()
274 { 283 {
275 // Header routine to compute various propert 284 // Header routine to compute various properties of material.
276 // << 285 //
>> 286
277 // Number of atoms per volume (per element), 287 // Number of atoms per volume (per element), total nb of electrons per volume
278 G4double Zi, Ai; 288 G4double Zi, Ai;
279 fTotNbOfAtomsPerVolume = 0.; << 289 TotNbOfAtomsPerVolume = 0.;
280 delete[] fVecNbOfAtomsPerVolume; << 290 if (VecNbOfAtomsPerVolume) { delete [] VecNbOfAtomsPerVolume; }
281 fVecNbOfAtomsPerVolume = new G4double[fNumbe << 291 VecNbOfAtomsPerVolume = new G4double[fNumberOfElements];
282 fTotNbOfElectPerVolume = 0.; << 292 TotNbOfElectPerVolume = 0.;
283 fFreeElecDensity = 0.; << 293 for (size_t i=0; i<fNumberOfElements; ++i) {
284 const G4double elecTh = 15. * CLHEP::eV; // << 294 Zi = (*theElementVector)[i]->GetZ();
285 for (G4int i = 0; i < fNumberOfElements; ++i << 295 Ai = (*theElementVector)[i]->GetA();
286 Zi = (*theElementVector)[i]->GetZ(); << 296 VecNbOfAtomsPerVolume[i] = Avogadro*fDensity*fMassFractionVector[i]/Ai;
287 Ai = (*theElementVector)[i]->GetA(); << 297 TotNbOfAtomsPerVolume += VecNbOfAtomsPerVolume[i];
288 fVecNbOfAtomsPerVolume[i] = Avogadro * fDe << 298 TotNbOfElectPerVolume += VecNbOfAtomsPerVolume[i]*Zi;
289 fTotNbOfAtomsPerVolume += fVecNbOfAtomsPer <<
290 fTotNbOfElectPerVolume += fVecNbOfAtomsPer <<
291 if (fState != kStateGas) { <<
292 fFreeElecDensity += <<
293 fVecNbOfAtomsPerVolume[i] * G4AtomicSh <<
294 } <<
295 } 299 }
296 << 300
297 ComputeRadiationLength(); 301 ComputeRadiationLength();
298 ComputeNuclearInterLength(); 302 ComputeNuclearInterLength();
299 303
300 if (fIonisation == nullptr) { << 304 if (fIonisation) { delete fIonisation; }
301 fIonisation = new G4IonisParamMat(this); << 305 fIonisation = new G4IonisParamMat(this);
302 } << 306 if (fSandiaTable) { delete fSandiaTable; }
303 if (fSandiaTable == nullptr) { << 307 fSandiaTable = new G4SandiaTable(this);
304 fSandiaTable = new G4SandiaTable(this); <<
305 } <<
306 } 308 }
307 309
308 //....oooOO0OOooo........oooOO0OOooo........oo 310 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
309 311
310 void G4Material::CopyPointersOfBaseMaterial() 312 void G4Material::CopyPointersOfBaseMaterial()
311 { 313 {
312 G4double factor = fDensity / fBaseMaterial-> << 314 G4double factor = fDensity/fBaseMaterial->GetDensity();
313 fTotNbOfAtomsPerVolume = factor * fBaseMater << 315 TotNbOfAtomsPerVolume = factor*fBaseMaterial->GetTotNbOfAtomsPerVolume();
314 fTotNbOfElectPerVolume = factor * fBaseMater << 316 TotNbOfElectPerVolume = factor*fBaseMaterial->GetTotNbOfElectPerVolume();
315 fFreeElecDensity = factor * fBaseMaterial->G <<
316 <<
317 if (fState == kStateUndefined) { <<
318 fState = fBaseMaterial->GetState(); <<
319 } <<
320 317
321 theElementVector = const_cast<G4ElementVecto 318 theElementVector = const_cast<G4ElementVector*>(fBaseMaterial->GetElementVector());
322 fMassFractionVector = const_cast<G4double*>( 319 fMassFractionVector = const_cast<G4double*>(fBaseMaterial->GetFractionVector());
323 fAtomsVector = const_cast<G4int*>(fBaseMater 320 fAtomsVector = const_cast<G4int*>(fBaseMaterial->GetAtomsVector());
324 321
325 const G4double* v = fBaseMaterial->GetVecNbO 322 const G4double* v = fBaseMaterial->GetVecNbOfAtomsPerVolume();
326 delete[] fVecNbOfAtomsPerVolume; << 323 if (VecNbOfAtomsPerVolume) { delete [] VecNbOfAtomsPerVolume; }
327 fVecNbOfAtomsPerVolume = new G4double[fNumbe << 324 VecNbOfAtomsPerVolume = new G4double[fNumberOfElements];
328 for (G4int i = 0; i < fNumberOfElements; ++i << 325 for (size_t i=0; i<fNumberOfElements; ++i) {
329 fVecNbOfAtomsPerVolume[i] = factor * v[i]; << 326 VecNbOfAtomsPerVolume[i] = factor*v[i];
330 } << 327 }
331 fRadlen = fBaseMaterial->GetRadlen() / facto << 328 fRadlen = fBaseMaterial->GetRadlen()/factor;
332 fNuclInterLen = fBaseMaterial->GetNuclearInt << 329 fNuclInterLen = fBaseMaterial->GetNuclearInterLength()/factor;
333 << 330 if (fIonisation) { delete fIonisation; }
334 if (fIonisation == nullptr) { << 331 fIonisation = new G4IonisParamMat(this);
335 fIonisation = new G4IonisParamMat(this); <<
336 } <<
337 fIonisation->SetMeanExcitationEnergy(fBaseMa <<
338 if (fBaseMaterial->GetIonisation()->GetDensi <<
339 ComputeDensityEffectOnFly(true); <<
340 } <<
341 332
342 fSandiaTable = fBaseMaterial->GetSandiaTable 333 fSandiaTable = fBaseMaterial->GetSandiaTable();
343 fMaterialPropertiesTable = fBaseMaterial->Ge << 334 fIonisation->SetMeanExcitationEnergy(fBaseMaterial->GetIonisation()->GetMeanExcitationEnergy());
344 } 335 }
345 336
346 //....oooOO0OOooo........oooOO0OOooo........oo 337 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
347 338
348 void G4Material::AddElementByNumberOfAtoms(con << 339 // AddElement -- composition by atom count
349 { <<
350 // perform checks consistency <<
351 if (0 == fIdxComponent) { <<
352 fMassFraction = false; <<
353 fAtoms = new std::vector<G4int>; <<
354 fElm = new std::vector<const G4Element*>; <<
355 } <<
356 if (fIdxComponent >= fNbComponents) { <<
357 G4ExceptionDescription ed; <<
358 ed << "For material " << fName << " and ad <<
359 << " with Natoms=" << nAtoms <<
360 << " wrong attempt to add more than the <<
361 << " >= " << fNbComponents; <<
362 G4Exception("G4Material::AddElementByNumbe <<
363 } <<
364 if (fMassFraction) { <<
365 G4ExceptionDescription ed; <<
366 ed << "For material " << fName << " and ad <<
367 << " with Natoms=" << nAtoms << " probl <<
368 << "addition of elements by mass fracti <<
369 G4Exception("G4Material::AddElementByNumbe <<
370 } <<
371 if (0 >= nAtoms) { <<
372 G4ExceptionDescription ed; <<
373 ed << "For material " << fName << " and ad <<
374 << " with Natoms=" << nAtoms << " probl <<
375 G4Exception("G4Material::AddElementByNumbe <<
376 } <<
377 <<
378 // filling <<
379 G4bool isAdded = false; <<
380 if (! fElm->empty()) { <<
381 for (G4int i = 0; i < fNumberOfElements; + <<
382 if (elm == (*fElm)[i]) { <<
383 (*fAtoms)[i] += nAtoms; <<
384 isAdded = true; <<
385 break; <<
386 } <<
387 } <<
388 } <<
389 if (! isAdded) { <<
390 fElm->push_back(elm); <<
391 fAtoms->push_back(nAtoms); <<
392 ++fNumberOfElements; <<
393 } <<
394 ++fIdxComponent; <<
395 <<
396 // is filled - complete composition of atoms <<
397 if (fIdxComponent == fNbComponents) { <<
398 theElementVector = new G4ElementVector(); <<
399 theElementVector->reserve(fNumberOfElement <<
400 fAtomsVector = new G4int[fNumberOfElements <<
401 fMassFractionVector = new G4double[fNumber <<
402 340
>> 341 void G4Material::AddElement(G4Element* element, G4int nAtoms)
>> 342 {
>> 343 // initialization
>> 344 if ( fNumberOfElements == 0 ) {
>> 345 fAtomsVector = new G4int [fArrayLength];
>> 346 fMassFractionVector = new G4double[fArrayLength];
>> 347 }
>> 348
>> 349 // filling ...
>> 350 if ( G4int(fNumberOfElements) < maxNbComponents ) {
>> 351 theElementVector->push_back(element);
>> 352 fAtomsVector[fNumberOfElements] = nAtoms;
>> 353 fNumberOfComponents = ++fNumberOfElements;
>> 354 element->increaseCountUse();
>> 355 } else {
>> 356 G4cout << "G4Material::AddElement ERROR for " << fName << " nElement= "
>> 357 << fNumberOfElements << G4endl;
>> 358 G4Exception ("G4Material::AddElement()", "mat031", FatalException,
>> 359 "Attempt to add more than the declared number of elements.");
>> 360 }
>> 361 // filled.
>> 362 if ( G4int(fNumberOfElements) == maxNbComponents ) {
>> 363 // compute proportion by mass
>> 364 size_t i=0;
403 G4double Amol = 0.; 365 G4double Amol = 0.;
404 for (G4int i = 0; i < fNumberOfElements; + << 366 for (i=0; i<fNumberOfElements; ++i) {
405 theElementVector->push_back((*fElm)[i]); << 367 G4double w = fAtomsVector[i]*(*theElementVector)[i]->GetA();
406 fAtomsVector[i] = (*fAtoms)[i]; <<
407 G4double w = fAtomsVector[i] * (*fElm)[i <<
408 Amol += w; 368 Amol += w;
409 fMassFractionVector[i] = w; 369 fMassFractionVector[i] = w;
410 } 370 }
411 for (G4int i = 0; i < fNumberOfElements; + << 371 for (i=0; i<fNumberOfElements; ++i) {
412 fMassFractionVector[i] /= Amol; 372 fMassFractionVector[i] /= Amol;
413 } 373 }
414 delete fAtoms; << 374
415 delete fElm; << 375 fMassOfMolecule = Amol/Avogadro;
416 fMassOfMolecule = Amol / CLHEP::Avogadro; <<
417 ComputeDerivedQuantities(); 376 ComputeDerivedQuantities();
418 } 377 }
419 } 378 }
420 379
421 //....oooOO0OOooo........oooOO0OOooo........oo 380 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
422 381
423 void G4Material::AddElementByMassFraction(cons << 382 // AddElement -- composition by fraction of mass
>> 383
>> 384 void G4Material::AddElement(G4Element* element, G4double fraction)
424 { 385 {
425 // perform checks consistency << 386 if(fraction < 0.0 || fraction > 1.0) {
426 if (fraction < 0.0 || fraction > 1.0) { << 387 G4cout << "G4Material::AddElement ERROR for " << fName << " and "
427 G4ExceptionDescription ed; << 388 << element->GetName() << " mass fraction= " << fraction
428 ed << "For material " << fName << " and ad << 389 << " is wrong " << G4endl;
429 << " massFraction= " << fraction << " i << 390 G4Exception ("G4Material::AddElement()", "mat032", FatalException,
430 G4Exception("G4Material::AddElementByMassF << 391 "Attempt to add element with wrong mass fraction");
431 } << 392 }
432 if (! fMassFraction) { << 393 // initialization
433 G4ExceptionDescription ed; << 394 if (fNumberOfComponents == 0) {
434 ed << "For material " << fName << " and ad << 395 fMassFractionVector = new G4double[fArrayLength];
435 << ", massFraction= " << fraction << ", << 396 fAtomsVector = new G4int [fArrayLength];
436 << " problem: cannot add by mass fracti << 397 }
437 << "addition of elements by number of a << 398 // filling ...
438 G4Exception("G4Material::AddElementByMassF << 399 if (G4int(fNumberOfComponents) < maxNbComponents) {
439 } << 400 size_t el = 0;
440 if (fIdxComponent >= fNbComponents) { << 401 while ((el<fNumberOfElements)&&(element!=(*theElementVector)[el])) { ++el; }
441 G4ExceptionDescription ed; << 402 if (el<fNumberOfElements) fMassFractionVector[el] += fraction;
442 ed << "For material " << fName << " and ad << 403 else {
443 << ", massFraction= " << fraction << ", << 404 theElementVector->push_back(element);
444 << "; attempt to add more than the decl << 405 fMassFractionVector[el] = fraction;
445 << " >= " << fNbComponents; << 406 ++fNumberOfElements;
446 G4Exception("G4Material::AddElementByMassF << 407 element->increaseCountUse();
447 } <<
448 if (0 == fIdxComponent) { <<
449 fElmFrac = new std::vector<G4double>; <<
450 fElm = new std::vector<const G4Element*>; <<
451 } <<
452 <<
453 // filling <<
454 G4bool isAdded = false; <<
455 if (! fElm->empty()) { <<
456 for (G4int i = 0; i < fNumberOfElements; + <<
457 if (elm == (*fElm)[i]) { <<
458 (*fElmFrac)[i] += fraction; <<
459 isAdded = true; <<
460 break; <<
461 } <<
462 } 408 }
463 } << 409 ++fNumberOfComponents;
464 if (! isAdded) { << 410 } else {
465 fElm->push_back(elm); << 411 G4cout << "G4Material::AddElement ERROR for " << fName << " nElement= "
466 fElmFrac->push_back(fraction); << 412 << fNumberOfElements << G4endl;
467 ++fNumberOfElements; << 413 G4Exception ("G4Material::AddElement()", "mat033", FatalException,
468 } << 414 "Attempt to add more than the declared number of elements.");
469 ++fIdxComponent; << 415 }
470 << 416
471 // is filled << 417 // filled.
472 if (fIdxComponent == fNbComponents) { << 418 if (G4int(fNumberOfComponents) == maxNbComponents) {
473 FillVectors(); << 419
>> 420 size_t i=0;
>> 421 G4double Zmol(0.), Amol(0.);
>> 422 // check sum of weights -- OK?
>> 423 G4double wtSum(0.0);
>> 424 for (i=0; i<fNumberOfElements; ++i) {
>> 425 wtSum += fMassFractionVector[i];
>> 426 Zmol += fMassFractionVector[i]*(*theElementVector)[i]->GetZ();
>> 427 Amol += fMassFractionVector[i]*(*theElementVector)[i]->GetA();
>> 428 }
>> 429 if (std::fabs(1.-wtSum) > perThousand) {
>> 430 G4cerr << "WARNING !! for " << fName << " sum of fractional masses "
>> 431 << wtSum << " is not 1 - results may be wrong"
>> 432 << G4endl;
>> 433 }
>> 434 for (i=0; i<fNumberOfElements; ++i) {
>> 435 fAtomsVector[i] =
>> 436 G4int(fMassFractionVector[i]*Amol/(*theElementVector)[i]->GetA()+0.5);
>> 437 }
>> 438
>> 439 ComputeDerivedQuantities();
474 } 440 }
475 } 441 }
476 442
477 //....oooOO0OOooo........oooOO0OOooo........oo 443 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
478 444
479 // composition by fraction of mass << 445 // AddMaterial -- composition by fraction of mass
>> 446
480 void G4Material::AddMaterial(G4Material* mater 447 void G4Material::AddMaterial(G4Material* material, G4double fraction)
481 { 448 {
482 if (fraction < 0.0 || fraction > 1.0) { << 449 if(fraction < 0.0 || fraction > 1.0) {
483 G4ExceptionDescription ed; << 450 G4cout << "G4Material::AddMaterial ERROR for " << fName << " and "
484 ed << "For material " << fName << " and ad << 451 << material->GetName() << " mass fraction= " << fraction
485 << ", massFraction= " << fraction << " << 452 << " is wrong " << G4endl;
486 G4Exception("G4Material::AddMaterial()", " << 453 G4Exception ("G4Material::AddMaterial()", "mat034", FatalException,
487 } << 454 "Attempt to add material with wrong mass fraction");
488 if (! fMassFraction) { << 455 }
489 G4ExceptionDescription ed; << 456 // initialization
490 ed << "For material " << fName << " and ad << 457 if (fNumberOfComponents == 0) {
491 << ", massFraction= " << fraction << ", << 458 fMassFractionVector = new G4double[fArrayLength];
492 << " problem: cannot add by mass fracti << 459 fAtomsVector = new G4int [fArrayLength];
493 << "addition of elements by number of a << 460 }
494 G4Exception("G4Material::AddMaterial()", " << 461
495 } << 462 size_t nelm = material->GetNumberOfElements();
496 if (fIdxComponent >= fNbComponents) { << 463
497 G4ExceptionDescription ed; << 464 // arrays should be extended
498 ed << "For material " << fName << " and ad << 465 if(nelm > 1) {
499 << ", massFraction= " << fraction << 466 G4int nold = fArrayLength;
500 << "; attempt to add more than the decl << 467 fArrayLength += nelm - 1;
501 << " >= " << fNbComponents; << 468 G4double* v1 = new G4double[fArrayLength];
502 G4Exception("G4Material::AddMaterial()", " << 469 G4int* i1 = new G4int[fArrayLength];
503 } << 470 for(G4int i=0; i<nold; ++i) {
504 if (0 == fIdxComponent) { << 471 v1[i] = fMassFractionVector[i];
505 fElmFrac = new std::vector<G4double>; << 472 i1[i] = fAtomsVector[i];
506 fElm = new std::vector<const G4Element*>; << 473 }
507 } << 474 delete [] fAtomsVector;
508 << 475 delete [] fMassFractionVector;
509 // filling << 476 fMassFractionVector = v1;
510 auto nelm = (G4int)material->GetNumberOfElem << 477 fAtomsVector = i1;
511 for (G4int j = 0; j < nelm; ++j) { << 478 }
512 auto elm = material->GetElement(j); << 479
513 auto frac = material->GetFractionVector(); << 480 // filling ...
514 G4bool isAdded = false; << 481 if (G4int(fNumberOfComponents) < maxNbComponents) {
515 if (! fElm->empty()) { << 482 for (size_t elm=0; elm<nelm; ++elm)
516 for (G4int i = 0; i < fNumberOfElements; << 483 {
517 if (elm == (*fElm)[i]) { << 484 G4Element* element = (*(material->GetElementVector()))[elm];
518 (*fElmFrac)[i] += fraction * frac[j] << 485 size_t el = 0;
519 isAdded = true; << 486 while ((el<fNumberOfElements)&&(element!=(*theElementVector)[el])) el++;
520 break; << 487 if (el < fNumberOfElements) fMassFractionVector[el] += fraction
>> 488 *(material->GetFractionVector())[elm];
>> 489 else {
>> 490 theElementVector->push_back(element);
>> 491 fMassFractionVector[el] = fraction
>> 492 *(material->GetFractionVector())[elm];
>> 493 ++fNumberOfElements;
>> 494 element->increaseCountUse();
521 } 495 }
522 } << 496 }
>> 497 ++fNumberOfComponents;
>> 498 ///store massFraction of material component
>> 499 fMatComponents[material] = fraction;
>> 500
>> 501 } else {
>> 502 G4cout << "G4Material::AddMaterial ERROR for " << fName << " nElement= "
>> 503 << fNumberOfElements << G4endl;
>> 504 G4Exception ("G4Material::AddMaterial()", "mat035", FatalException,
>> 505 "Attempt to add more than the declared number of components.");
>> 506 }
>> 507
>> 508 // filled.
>> 509 if (G4int(fNumberOfComponents) == maxNbComponents) {
>> 510 size_t i=0;
>> 511 G4double Zmol(0.), Amol(0.);
>> 512 // check sum of weights -- OK?
>> 513 G4double wtSum(0.0);
>> 514 for (i=0; i<fNumberOfElements; ++i) {
>> 515 wtSum += fMassFractionVector[i];
>> 516 Zmol += fMassFractionVector[i]*(*theElementVector)[i]->GetZ();
>> 517 Amol += fMassFractionVector[i]*(*theElementVector)[i]->GetA();
>> 518 }
>> 519 if (std::fabs(1.-wtSum) > perThousand) {
>> 520 G4cout << "G4Material::AddMaterial WARNING !! for " << fName
>> 521 << " sum of fractional masses "
>> 522 << wtSum << " is not 1 - results may be wrong"
>> 523 << G4endl;
>> 524 }
>> 525 for (i=0;i<fNumberOfElements;i++) {
>> 526 fAtomsVector[i] =
>> 527 G4int(fMassFractionVector[i]*Amol/(*theElementVector)[i]->GetA()+0.5);
523 } 528 }
524 if (! isAdded) { << 529
525 fElm->push_back(elm); << 530 ComputeDerivedQuantities();
526 fElmFrac->push_back(fraction * frac[j]); <<
527 ++fNumberOfElements; <<
528 } <<
529 } <<
530 <<
531 fMatComponents[material] = fraction; <<
532 ++fIdxComponent; <<
533 <<
534 // is filled <<
535 if (fIdxComponent == fNbComponents) { <<
536 FillVectors(); <<
537 } <<
538 } <<
539 <<
540 //....oooOO0OOooo........oooOO0OOooo........oo <<
541 <<
542 void G4Material::FillVectors() <<
543 { <<
544 // there are material components <<
545 theElementVector = new G4ElementVector(); <<
546 theElementVector->reserve(fNumberOfElements) <<
547 fAtomsVector = new G4int[fNumberOfElements]; <<
548 fMassFractionVector = new G4double[fNumberOf <<
549 <<
550 G4double wtSum(0.0); <<
551 for (G4int i = 0; i < fNumberOfElements; ++i <<
552 theElementVector->push_back((*fElm)[i]); <<
553 fMassFractionVector[i] = (*fElmFrac)[i]; <<
554 wtSum += fMassFractionVector[i]; <<
555 } <<
556 delete fElmFrac; <<
557 delete fElm; <<
558 <<
559 // check sum of weights -- OK? <<
560 if (std::abs(1. - wtSum) > perThousand) { <<
561 G4ExceptionDescription ed; <<
562 ed << "For material " << fName << " sum of <<
563 << " is not 1 - results may be wrong"; <<
564 G4Exception("G4Material::FillVectors()", " <<
565 } <<
566 G4double coeff = (wtSum > 0.0) ? 1. / wtSum <<
567 G4double Amol(0.); <<
568 for (G4int i = 0; i < fNumberOfElements; ++i <<
569 fMassFractionVector[i] *= coeff; <<
570 Amol += fMassFractionVector[i] * (*theElem <<
571 } <<
572 for (G4int i = 0; i < fNumberOfElements; ++i <<
573 fAtomsVector[i] = G4lrint(fMassFractionVec <<
574 } 531 }
575 ComputeDerivedQuantities(); <<
576 } 532 }
577 533
578 //....oooOO0OOooo........oooOO0OOooo........oo 534 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
579 535
580 void G4Material::ComputeRadiationLength() 536 void G4Material::ComputeRadiationLength()
581 { 537 {
582 G4double radinv = 0.0; << 538 G4double radinv = 0.0 ;
583 for (G4int i = 0; i < fNumberOfElements; ++i << 539 for (size_t i=0;i<fNumberOfElements;++i) {
584 radinv += fVecNbOfAtomsPerVolume[i] * ((*t << 540 radinv += VecNbOfAtomsPerVolume[i]*((*theElementVector)[i]->GetfRadTsai());
585 } 541 }
586 fRadlen = (radinv <= 0.0 ? DBL_MAX : 1. / ra << 542 fRadlen = (radinv <= 0.0 ? DBL_MAX : 1./radinv);
587 } 543 }
588 544
589 //....oooOO0OOooo........oooOO0OOooo........oo 545 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
590 546
591 void G4Material::ComputeNuclearInterLength() 547 void G4Material::ComputeNuclearInterLength()
592 { 548 {
593 const G4double lambda1 = CLHEP::amu*CLHEP::c << 549 const G4double lambda0 = 35*g/cm2;
594 G4double NILinv = 0.0; 550 G4double NILinv = 0.0;
595 for (G4int i = 0; i < fNumberOfElements; ++i << 551 G4Pow* g4pow = G4Pow::GetInstance();
596 G4int Z = (*theElementVector)[i]->GetZasIn << 552 for (size_t i=0; i<fNumberOfElements; ++i) {
597 G4double A = (*theElementVector)[i]->GetN( << 553 NILinv +=
598 if (1 == Z) { << 554 VecNbOfAtomsPerVolume[i]*g4pow->Z23(G4int((*theElementVector)[i]->GetN()+0.5));
599 NILinv += fVecNbOfAtomsPerVolume[i] * A; <<
600 } <<
601 else { <<
602 NILinv += fVecNbOfAtomsPerVolume[i] * G4 <<
603 } <<
604 } <<
605 NILinv *= lambda1; <<
606 fNuclInterLen = 1.e+20*CLHEP::m; <<
607 if (fNuclInterLen * NILinv > 1.0) { <<
608 fNuclInterLen = 1.0 / NILinv; <<
609 } 555 }
>> 556 NILinv *= amu/lambda0;
>> 557 fNuclInterLen = (NILinv <= 0.0 ? DBL_MAX : 1./NILinv);
610 } 558 }
611 559
612 //....oooOO0OOooo........oooOO0OOooo........oo 560 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
613 561
614 void G4Material::SetChemicalFormula(const G4St << 562 const G4MaterialTable* G4Material::GetMaterialTable()
615 { 563 {
616 if (! IsLocked()) { << 564 return &theMaterialTable;
617 fChemicalFormula = chF; <<
618 } <<
619 } 565 }
620 566
621 //....oooOO0OOooo........oooOO0OOooo........oo 567 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
622 568
623 void G4Material::SetFreeElectronDensity(G4doub << 569 size_t G4Material::GetNumberOfMaterials()
624 { 570 {
625 if (val >= 0. && ! IsLocked()) { << 571 return theMaterialTable.size();
626 fFreeElecDensity = val; <<
627 } <<
628 } 572 }
629 573
630 //....oooOO0OOooo........oooOO0OOooo........oo 574 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
631 575
632 void G4Material::ComputeDensityEffectOnFly(G4b << 576 G4Material* G4Material::GetMaterial(const G4String& materialName, G4bool warning)
633 { << 577 {
634 if (! IsLocked()) { << 578 // search the material by its name
635 if (nullptr == fIonisation) { << 579 for (size_t J=0 ; J<theMaterialTable.size() ; ++J)
636 fIonisation = new G4IonisParamMat(this); << 580 {
>> 581 if (theMaterialTable[J]->GetName() == materialName)
>> 582 { return theMaterialTable[J]; }
637 } 583 }
638 fIonisation->ComputeDensityEffectOnFly(val << 584
639 } << 585 // the material does not exist in the table
>> 586 if (warning) {
>> 587 G4cout << "G4Material::GetMaterial() WARNING: The material: "
>> 588 << materialName << " does not exist in the table. Return NULL pointer."
>> 589 << G4endl;
>> 590 }
>> 591 return 0;
640 } 592 }
641 593
642 //....oooOO0OOooo........oooOO0OOooo........oo 594 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
643 595
644 G4MaterialTable* G4Material::GetMaterialTable( << 596 G4Material::G4Material(const G4Material& right)
645 { 597 {
646 struct Holder { << 598 InitializePointers();
647 G4MaterialTable instance; << 599 *this = right;
648 ~Holder() { <<
649 for(auto item : instance) <<
650 delete item; <<
651 } <<
652 }; <<
653 static Holder _holder; <<
654 return &_holder.instance; <<
655 } 600 }
656 601
657 //....oooOO0OOooo........oooOO0OOooo........oo 602 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
658 603
659 std::size_t G4Material::GetNumberOfMaterials() << 604 G4double G4Material::GetZ() const
>> 605 {
>> 606 if (fNumberOfElements > 1) {
>> 607 G4cout << "G4Material ERROR in GetZ. The material: " << fName << " is a mixture."
>> 608 << G4endl;
>> 609 G4Exception ("G4Material::GetZ()", "mat036", FatalException,
>> 610 "the Atomic number is not well defined." );
>> 611 }
>> 612 return (*theElementVector)[0]->GetZ();
>> 613 }
660 614
661 //....oooOO0OOooo........oooOO0OOooo........oo 615 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
662 616
663 G4Material* G4Material::GetMaterial(const G4St << 617 G4double G4Material::GetA() const
664 { << 618 {
665 // search the material by its name << 619 if (fNumberOfElements > 1) {
666 for (auto const & j : *GetMaterialTable()) { << 620 G4cout << "G4Material ERROR in GetA. The material: " << fName << " is a mixture."
667 if (j->GetName() == materialName) { << 621 << G4endl;
668 return j; << 622 G4Exception ("G4Material::GetA()", "mat037", FatalException,
669 } << 623 "the Atomic mass is not well defined." );
670 } << 624 }
671 << 625 return (*theElementVector)[0]->GetA();
672 // the material does not exist in the table <<
673 if (warn) { <<
674 G4cout << "G4Material::GetMaterial() WARNI <<
675 << " does not exist in the table. R <<
676 } <<
677 return nullptr; <<
678 } 626 }
679 627
680 //....oooOO0OOooo........oooOO0OOooo........oo 628 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
681 629
682 G4Material* G4Material::GetMaterial(G4double z << 630 const G4Material& G4Material::operator=(const G4Material& right)
683 { 631 {
684 // search the material by its name << 632 if (this != &right)
685 for (auto const & mat : *GetMaterialTable()) <<
686 if (1 == mat->GetNumberOfElements() && z = <<
687 dens == mat->GetDensity()) <<
688 { 633 {
689 return mat; << 634 fName = right.fName;
690 } << 635 fChemicalFormula = right.fChemicalFormula;
691 } << 636 fDensity = right.fDensity;
692 return nullptr; << 637 fState = right.fState;
693 } << 638 fTemp = right.fTemp;
>> 639 fPressure = right.fPressure;
>> 640
>> 641 if(!fBaseMaterial) {
>> 642 if (theElementVector) { delete theElementVector; }
>> 643 if (fMassFractionVector) { delete [] fMassFractionVector; }
>> 644 if (fAtomsVector) { delete [] fAtomsVector; }
>> 645 if (fIonisation) { delete fIonisation; }
>> 646 if (fSandiaTable) { delete fSandiaTable; }
>> 647 }
694 648
695 //....oooOO0OOooo........oooOO0OOooo........oo << 649 if (VecNbOfAtomsPerVolume) { delete [] VecNbOfAtomsPerVolume; }
696 650
697 G4Material* G4Material::GetMaterial(std::size_ << 651 maxNbComponents = right.maxNbComponents;
698 { << 652 fNumberOfComponents = right.fNumberOfComponents;
699 // search the material by its name << 653 fNumberOfElements = right.fNumberOfElements;
700 for (auto const & mat : *GetMaterialTable()) << 654 fImplicitElement = right.fImplicitElement;
701 if (nComp == mat->GetNumberOfElements() && << 655
702 return mat; << 656 fMaterialPropertiesTable = right.fMaterialPropertiesTable;
703 } << 657 fBaseMaterial = right.fBaseMaterial;
704 } << 658 fMassOfMolecule= right.fMassOfMolecule;
705 return nullptr; << 659 fMatComponents= right.fMatComponents;
>> 660
>> 661 if(fBaseMaterial) {
>> 662 CopyPointersOfBaseMaterial();
>> 663
>> 664 } else {
>> 665 theElementVector = new G4ElementVector(fNumberOfElements,0);
>> 666 fMassFractionVector = new G4double[fNumberOfElements];
>> 667 fAtomsVector = new G4int[fNumberOfElements];
>> 668 for (size_t i=0; i<fNumberOfElements; ++i) {
>> 669 (*theElementVector)[i] = (*right.theElementVector)[i];
>> 670 fMassFractionVector[i] = right.fMassFractionVector[i];
>> 671 fAtomsVector[i] = right.fAtomsVector[i];
>> 672 }
>> 673 ComputeDerivedQuantities();
>> 674 }
>> 675
>> 676 }
>> 677 return *this;
706 } 678 }
707 679
708 //....oooOO0OOooo........oooOO0OOooo........oo 680 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
709 681
710 G4double G4Material::GetZ() const << 682 G4int G4Material::operator==(const G4Material& right) const
711 { 683 {
712 if (fNumberOfElements > 1) { << 684 return (this == (G4Material *) &right);
713 G4ExceptionDescription ed; <<
714 ed << "For material " << fName << " ERROR <<
715 << " > 1, which is not allowed"; <<
716 G4Exception("G4Material::GetZ()", "mat036" <<
717 } <<
718 return (*theElementVector)[0]->GetZ(); <<
719 } 685 }
720 686
721 //....oooOO0OOooo........oooOO0OOooo........oo 687 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
722 688
723 G4double G4Material::GetA() const << 689 G4int G4Material::operator!=(const G4Material& right) const
724 { 690 {
725 if (fNumberOfElements > 1) { << 691 return (this != (G4Material *) &right);
726 G4ExceptionDescription ed; <<
727 ed << "For material " << fName << " ERROR <<
728 << " > 1, which is not allowed"; <<
729 G4Exception("G4Material::GetA()", "mat036" <<
730 } <<
731 return (*theElementVector)[0]->GetA(); <<
732 } 692 }
733 693
734 //....oooOO0OOooo........oooOO0OOooo........oo 694 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
735 695
736 std::ostream& operator<<(std::ostream& flux, c << 696
>> 697 std::ostream& operator<<(std::ostream& flux, G4Material* material)
737 { 698 {
738 std::ios::fmtflags mode = flux.flags(); 699 std::ios::fmtflags mode = flux.flags();
739 flux.setf(std::ios::fixed, std::ios::floatfi << 700 flux.setf(std::ios::fixed,std::ios::floatfield);
740 G4long prec = flux.precision(3); 701 G4long prec = flux.precision(3);
741 << 702
742 flux << " Material: " << std::setw(8) << mat << 703 flux
743 << " " << 704 << " Material: " << std::setw(8) << material->fName
744 << " density: " << std::setw(6) << std << 705 << " " << material->fChemicalFormula << " "
745 << G4BestUnit(material->fDensity, "Volu << 706 << " density: " << std::setw(6) << std::setprecision(3)
746 << std::setprecision(3) << G4BestUnit(m << 707 << G4BestUnit(material->fDensity,"Volumic Mass")
747 << " Nucl.Int.Length: " << std::setw(7 << 708 << " RadL: " << std::setw(7) << std::setprecision(3)
748 << G4BestUnit(material->fNuclInterLen, << 709 << G4BestUnit(material->fRadlen,"Length")
749 << std::setw(30) << " Imean: " << std: << 710 << " Nucl.Int.Length: " << std::setw(7) << std::setprecision(3)
750 << G4BestUnit(material->GetIonisation() << 711 << G4BestUnit(material->fNuclInterLen,"Length")
751 << " temperature: " << std::setw(6) << << 712 << " Imean: " << std::setw(7) << std::setprecision(3)
752 << (material->fTemp) / CLHEP::kelvin << << 713 << G4BestUnit(material->GetIonisation()->GetMeanExcitationEnergy(),"Energy");
753 << " pressure: " << std::setw(6) << st << 714
754 << (material->fPressure) / CLHEP::atmos << 715 if(material->fState == kStateGas) {
755 << "\n"; << 716 flux
756 << 717 << " temperature: " << std::setw(6) << std::setprecision(2)
757 for (G4int i = 0; i < material->fNumberOfEle << 718 << (material->fTemp)/kelvin << " K"
758 flux << "\n ---> " << (*(material->theEl << 719 << " pressure: " << std::setw(6) << std::setprecision(2)
759 << "\n ElmMassFraction: " << << 720 << (material->fPressure)/atmosphere << " atm";
760 << (material->fMassFractionVector[i]) << 721 }
761 << " ElmAbundance " << std::setw(6) << 722 for (size_t i=0; i<material->fNumberOfElements; i++) {
762 << 100 * (material->fVecNbOfAtomsPerV << 723 flux
763 << " % \n"; << 724 << "\n ---> " << (*(material->theElementVector))[i]
764 } << 725 << "\n ElmMassFraction: "
765 flux.precision(prec); << 726 << std::setw(6)<< std::setprecision(2)
766 flux.setf(mode, std::ios::floatfield); << 727 << (material->fMassFractionVector[i])/perCent << " %"
767 << 728 << " ElmAbundance " << std::setw(6)<< std::setprecision(2)
768 if (material->IsExtended()) { << 729 << 100*(material->VecNbOfAtomsPerVolume[i])/(material->TotNbOfAtomsPerVolume)
769 static_cast<const G4ExtendedMaterial*>(mat << 730 << " % \n";
770 } << 731 }
771 << 732 flux.precision(prec);
>> 733 flux.setf(mode,std::ios::floatfield);
>> 734
772 return flux; 735 return flux;
773 } 736 }
774 737
775 //....oooOO0OOooo........oooOO0OOooo........oo 738 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
776 739
777 std::ostream& operator<<(std::ostream& flux, c << 740 std::ostream& operator<<(std::ostream& flux, G4Material& material)
778 { 741 {
779 flux << &material; << 742 flux << &material;
780 return flux; 743 return flux;
781 } 744 }
782 745
783 //....oooOO0OOooo........oooOO0OOooo........oo 746 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
784 << 747
785 std::ostream& operator<<(std::ostream& flux, c << 748 std::ostream& operator<<(std::ostream& flux, G4MaterialTable MaterialTable)
786 { 749 {
787 // Dump info for all known materials << 750 //Dump info for all known materials
788 flux << "\n***** Table : Nb of materials = " << 751 flux << "\n***** Table : Nb of materials = " << MaterialTable.size()
789 << 752 << " *****\n" << G4endl;
790 for (auto i : MaterialTable) { << 753
791 flux << i << G4endl << G4endl; << 754 for (size_t i=0; i<MaterialTable.size(); ++i) {
>> 755 flux << MaterialTable[i] << G4endl << G4endl;
792 } 756 }
793 757
794 return flux; 758 return flux;
795 } << 759 }
796 <<
797 //....oooOO0OOooo........oooOO0OOooo........oo <<
798 <<
799 G4bool G4Material::IsExtended() const { return <<
800 760
801 //....oooOO0OOooo........oooOO0OOooo........oo 761 //....oooOO0OOooo........oooOO0OOooo........oooOO0OOooo........oooOO0OOooo......
802 <<
803 void G4Material::SetMaterialPropertiesTable(G4 <<
804 { <<
805 if (fMaterialPropertiesTable != anMPT && ! I <<
806 delete fMaterialPropertiesTable; <<
807 fMaterialPropertiesTable = anMPT; <<
808 } <<
809 } <<
810 <<
811 //....oooOO0OOooo........oooOO0OOooo........oo <<
812 <<
813 G4bool G4Material::IsLocked() <<
814 { <<
815 auto state = G4StateManager::GetStateManager <<
816 return state != G4State_PreInit && state != <<
817 } <<
818 762