| Geant4 Cross Reference |
1 // 1
2 // *******************************************
3 // * License and Disclaimer
4 // *
5 // * The Geant4 software is copyright of th
6 // * the Geant4 Collaboration. It is provided
7 // * conditions of the Geant4 Software License
8 // * LICENSE and available at http://cern.ch/
9 // * include a list of copyright holders.
10 // *
11 // * Neither the authors of this software syst
12 // * institutes,nor the agencies providing fin
13 // * work make any representation or warran
14 // * regarding this software system or assum
15 // * use. Please see the license in the file
16 // * for the full disclaimer and the limitatio
17 // *
18 // * This code implementation is the result
19 // * technical work of the GEANT4 collaboratio
20 // * By using, copying, modifying or distri
21 // * any work based on the software) you ag
22 // * use in resulting scientific publicati
23 // * acceptance of all terms of the Geant4 Sof
24 // *******************************************
25 //
26 //
27 // -------------------------------------------
28 //
29 // Geant4 source file
30 //
31 // File name: G4ParticleHPIsoProbabilityT
32 //
33 // Authors: Marek Zmeskal (CTU, Czech Tec
34 // Loic Thulliez (CEA France)
35 //
36 // Creation date: 4 June 2024
37 //
38 // Description: Class for the probability
39 // and for the given tempera
40 // It reads the files with p
41 // finds the correct cross-s
42 //
43 // Modifications:
44 //
45 // -------------------------------------------
46 //
47 //
48
49 #include "G4ParticleHPIsoProbabilityTable.hh"
50 #include "G4SystemOfUnits.hh"
51 #include "Randomize.hh"
52 #include "G4Exception.hh"
53 #include "G4NucleiProperties.hh"
54 #include "G4ReactionProduct.hh"
55 #include "G4ParticleHPManager.hh"
56 #include "G4ParticleHPVector.hh"
57 #include "G4DynamicParticle.hh"
58 #include "G4Element.hh"
59 #include "G4Nucleus.hh"
60 #include "G4Neutron.hh"
61 #include "G4ParticleHPChannel.hh"
62 #include "G4ParticleHPChannelList.hh"
63
64 #include <string>
65 #include <fstream>
66
67 ///-------------------------------------------
68 G4ParticleHPIsoProbabilityTable::~G4ParticleHP
69 {
70 for ( auto it = theProbabilities->cbegin();
71 delete* it;
72 }
73 for ( auto it = theElasticData->cbegin(); it
74 delete* it;
75 }
76 for ( auto it = theCaptureData->cbegin(); it
77 delete* it;
78 }
79 for ( auto it = theFissionData->cbegin(); it
80 delete* it;
81 }
82 delete theEnergies;
83 delete theProbabilities;
84 delete theElasticData;
85 delete theCaptureData;
86 delete theFissionData;
87 }
88
89 ///-------------------------------------------
90 void G4ParticleHPIsoProbabilityTable::Init( G4
91
92 ///-------------------------------------------
93 G4double G4ParticleHPIsoProbabilityTable::GetC
94 const G4Element*, G4double&, G4double&, std:
95 {
96 G4Exception( "G4ParticleHPIsoProbabilityTabl
97 FatalException, "The base method G4Pa
98 "is trying to be accessed, which is n
99 "the derived class for NJOY or CALEND
100 return 0.0;
101 }
102
103 ///-------------------------------------------
104 G4double G4ParticleHPIsoProbabilityTable::GetI
105 const G4Element*, G4double&, std::map< std::
106 {
107 G4Exception( "G4ParticleHPIsoProbabilityTabl
108 FatalException, "The base method G4Pa
109 "is trying to be accessed, which is n
110 "the derived class for NJOY or CALEND
111 return 0.0;
112 }
113
114 ///-------------------------------------------
115 G4double G4ParticleHPIsoProbabilityTable::GetD
116 G4double result = 0.0;
117 G4ReactionProduct theNeutron( dP->GetDefinit
118 theNeutron.SetMomentum( dP->GetMomentum() );
119 theNeutron.SetKineticEnergy( dP->GetKineticE
120 // prepare thermal nucleus
121 G4Nucleus aNuc;
122 G4double eleMass = G4NucleiProperties::GetNu
123 G4ReactionProduct boosted;
124 G4double aXsection;
125 G4int counter = 0;
126 G4double buffer = 0.0;
127 G4int size = G4int( std::max( 10.0, T / 60.0
128 G4ThreeVector neutronVelocity = 1.0 / G4Neut
129 G4double neutronVMag = neutronVelocity.mag()
130 while ( counter == 0 || std::abs( buffer -
131 if ( counter ) buffer = result / counter;
132 while ( counter < size ) {
133 ++counter;
134 G4ReactionProduct aThermalNuc = aNuc.Get
135 boosted.Lorentz( theNeutron, aThermalNuc
136 G4double theEkin = boosted.GetKineticEne
137 aXsection = (*G4ParticleHPManager::GetIn
138 // velocity correction
139 G4ThreeVector targetVelocity = 1.0 / aTh
140 aXsection *= (targetVelocity - neutronVe
141 result += aXsection;
142 }
143 size += size;
144 }
145 result /= counter;
146 return result;
147 }
148
149 ///-------------------------------------------
150 G4double G4ParticleHPIsoProbabilityTable::GetD
151 G4double result = 0.0;
152 G4ReactionProduct theNeutron( dP->GetDefinit
153 theNeutron.SetMomentum( dP->GetMomentum() );
154 theNeutron.SetKineticEnergy( dP->GetKineticE
155 // prepare thermal nucleus
156 G4Nucleus aNuc;
157 G4double eleMass = G4NucleiProperties::GetNu
158 G4ReactionProduct boosted;
159 G4double aXsection;
160 G4int counter = 0;
161 G4double buffer = 0.0;
162 G4int size = G4int( std::max( 10.0, T / 60.0
163 G4ThreeVector neutronVelocity = 1.0 / G4Neut
164 G4double neutronVMag = neutronVelocity.mag()
165 while ( counter == 0 || std::abs( buffer -
166 if ( counter ) buffer = result / counter;
167 while ( counter < size ) {
168 ++counter;
169 G4ReactionProduct aThermalNuc = aNuc.Get
170 boosted.Lorentz( theNeutron, aThermalNuc
171 G4double theEkin = boosted.GetKineticEne
172 aXsection = (*G4ParticleHPManager::GetIn
173 // velocity correction
174 G4ThreeVector targetVelocity = 1.0 / aTh
175 aXsection *= (targetVelocity - neutronVe
176 result += aXsection;
177 }
178 size += size;
179 }
180 result /= counter;
181 return result;
182 }
183
184 ///-------------------------------------------
185 G4double G4ParticleHPIsoProbabilityTable::GetD
186 G4double result = 0.0;
187 G4ReactionProduct theNeutron( dP->GetDefinit
188 theNeutron.SetMomentum( dP->GetMomentum() );
189 theNeutron.SetKineticEnergy( dP->GetKineticE
190 // prepare thermal nucleus
191 G4Nucleus aNuc;
192 G4double eleMass;
193 eleMass = G4NucleiProperties::GetNuclearMass
194 G4ReactionProduct boosted;
195 G4double aXsection;
196 G4int counter = 0;
197 G4double buffer = 0.0;
198 G4int size = G4int( std::max( 10.0, T / 60.0
199 G4ThreeVector neutronVelocity = 1.0 / G4Neut
200 G4double neutronVMag = neutronVelocity.mag()
201 while ( counter == 0 || std::abs( buffer -
202 if ( counter ) buffer = result / counter;
203 while ( counter < size ) {
204 ++counter;
205 G4ReactionProduct aThermalNuc = aNuc.Get
206 boosted.Lorentz( theNeutron, aThermalNuc
207 G4double theEkin = boosted.GetKineticEne
208 aXsection = (*G4ParticleHPManager::GetIn
209 // velocity correction
210 G4ThreeVector targetVelocity = 1.0 / aTh
211 aXsection *= (targetVelocity - neutronVe
212 result += aXsection;
213 }
214 size += size;
215 }
216 result /= counter;
217 return result;
218 }
219
220 ///-------------------------------------------
221 G4double G4ParticleHPIsoProbabilityTable::GetD
222 G4double result = 0.0;
223 G4ReactionProduct theNeutron( dP->GetDefinit
224 theNeutron.SetMomentum( dP->GetMomentum() );
225 theNeutron.SetKineticEnergy( dP->GetKineticE
226 // prepare thermal nucleus
227 G4Nucleus aNuc;
228 G4double eleMass = G4NucleiProperties::GetNu
229 G4ReactionProduct boosted;
230 G4double aXsection;
231 G4int counter = 0;
232 G4int failCount = 0;
233 G4double buffer = 0.0;
234 G4int size = G4int( std::max( 10.0, T / 60.0
235 G4ThreeVector neutronVelocity = 1.0 / G4Neut
236 G4double neutronVMag = neutronVelocity.mag()
237 #ifndef G4PHP_DOPPLER_LOOP_ONCE
238 while ( counter == 0 || std::abs( buffer -
239 if ( counter ) buffer = result / counter;
240 while ( counter < size ) {
241 ++counter;
242 #endif
243 G4ReactionProduct aThermalNuc = aNuc.Get
244 boosted.Lorentz( theNeutron, aThermalNuc
245 G4double theEkin = boosted.GetKineticEne
246 aXsection = ((*G4ParticleHPManager::GetI
247 ->GetWeightedXsec( theEkin,
248 if ( aXsection < 0.0 ) {
249 if ( failCount < 1000 ) {
250 ++failCount;
251 #ifndef G4PHP_DOPPLER_LOOP_ONCE
252 --counter;
253 continue;
254 #endif
255 } else {
256 aXsection = 0.0;
257 }
258 }
259 // velocity correction
260 G4ThreeVector targetVelocity = 1.0 / aTh
261 aXsection *= (targetVelocity - neutronVe
262 result += aXsection;
263 }
264 #ifndef G4PHP_DOPPLER_LOOP_ONCE
265 size += size;
266 }
267 result /= counter;
268 #endif
269 return result;
270 }
271