| 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 // INCL++ intra-nuclear cascade model 26 // INCL++ intra-nuclear cascade model
27 // Alain Boudard, CEA-Saclay, France << 27 // Pekka Kaitaniemi, CEA and Helsinki Institute of Physics
28 // Joseph Cugnon, University of Liege, Belgium << 28 // Davide Mancusi, CEA
29 // Jean-Christophe David, CEA-Saclay, France << 29 // Alain Boudard, CEA
30 // Pekka Kaitaniemi, CEA-Saclay, France, and H << 30 // Sylvie Leray, CEA
31 // Sylvie Leray, CEA-Saclay, France << 31 // Joseph Cugnon, University of Liege
32 // Davide Mancusi, CEA-Saclay, France << 32 //
>> 33 // INCL++ revision: v5.0_rc3
33 // 34 //
34 #define INCLXX_IN_GEANT4_MODE 1 35 #define INCLXX_IN_GEANT4_MODE 1
35 36
36 #include "globals.hh" 37 #include "globals.hh"
37 38
38 #ifndef G4INCLClusteringModelIntercomparison_h 39 #ifndef G4INCLClusteringModelIntercomparison_hh
39 #define G4INCLClusteringModelIntercomparison_h 40 #define G4INCLClusteringModelIntercomparison_hh 1
40 41
41 #ifdef INCLXX_IN_GEANT4_MODE <<
42 #define INCL_CACHING_CLUSTERING_MODEL_INTERCOM <<
43 #endif // INCLXX_IN_GEANT4_MODE <<
44 <<
45 #include "G4INCLIClusteringModel.hh" 42 #include "G4INCLIClusteringModel.hh"
46 #include "G4INCLParticle.hh" 43 #include "G4INCLParticle.hh"
47 #include "G4INCLParticleTable.hh" <<
48 #include "G4INCLCluster.hh" 44 #include "G4INCLCluster.hh"
49 #include "G4INCLNucleus.hh" 45 #include "G4INCLNucleus.hh"
50 #include "G4INCLKinematicsUtils.hh" 46 #include "G4INCLKinematicsUtils.hh"
51 #include "G4INCLHashing.hh" <<
52 <<
53 #include <set> <<
54 #include <algorithm> <<
55 47
56 namespace G4INCL { 48 namespace G4INCL {
57 49
58 /** \brief Container for the relevant inform <<
59 * <<
60 * This struct contains all the information <<
61 * clustering algorithm. It is probably more <<
62 * feeds on, hopefully improving cache perfo <<
63 */ <<
64 struct ConsideredPartner { <<
65 Particle *particle; <<
66 G4bool isTargetSpectator; <<
67 G4int Z; <<
68 G4int S; <<
69 ThreeVector position; <<
70 ThreeVector momentum; <<
71 G4double energy; <<
72 G4double potentialEnergy; <<
73 <<
74 ConsideredPartner() : <<
75 particle(NULL), <<
76 isTargetSpectator(false), <<
77 Z(0), <<
78 S(0), <<
79 energy(0.), <<
80 potentialEnergy(0.) <<
81 {} <<
82 <<
83 ConsideredPartner(Particle * const p) : <<
84 particle(p), <<
85 isTargetSpectator(particle->isTargetSpec <<
86 Z(particle->getZ()), <<
87 S(particle->getS()), <<
88 position(particle->getPosition()), <<
89 momentum(particle->getMomentum()), <<
90 energy(particle->getEnergy()), <<
91 potentialEnergy(particle->getPotentialEn <<
92 {} <<
93 }; <<
94 <<
95 /// \brief Cluster coalescence algorithm use <<
96 class ClusteringModelIntercomparison : publi 50 class ClusteringModelIntercomparison : public IClusteringModel {
97 public: 51 public:
98 ClusteringModelIntercomparison(Config cons << 52 ClusteringModelIntercomparison() {
99 theNucleus(NULL), << 53 zeroOut();
100 selectedA(0), << 54
101 selectedZ(0), <<
102 selectedS(0), <<
103 sqtot(0.), <<
104 cascadingEnergyPool(0.), <<
105 protonMass(ParticleTable::getRealMass(Pr <<
106 neutronMass(ParticleTable::getRealMass(N <<
107 lambdaMass(ParticleTable::getRealMass(La <<
108 runningMaxClusterAlgorithmMass(theConfig <<
109 nConsideredMax(0), <<
110 nConsidered(0), <<
111 consideredPartners(NULL), <<
112 isInRunningConfiguration(NULL), <<
113 maxMassConfigurationSkipping(ParticleTab <<
114 { <<
115 // Set up the maximum charge and neutron 55 // Set up the maximum charge and neutron number for clusters
116 clusterZMaxAll = 0; 56 clusterZMaxAll = 0;
117 clusterNMaxAll = 0; 57 clusterNMaxAll = 0;
118 for(G4int A=0; A<=runningMaxClusterAlgor << 58 for(G4int A=0; A<=maxClusterAlgorithmMass; ++A) {
119 if(clusterZMax[A]>clusterZMaxAll) << 59 if(ParticleTable::clusterZMax[A]>clusterZMaxAll)
120 clusterZMaxAll = clusterZMax[A]; << 60 clusterZMaxAll = ParticleTable::clusterZMax[A];
121 if(A-clusterZMin[A]>clusterNMaxAll) << 61 if(A-ParticleTable::clusterZMin[A]>clusterNMaxAll)
122 clusterNMaxAll = A-clusterZMin[A]; << 62 clusterNMaxAll = A-ParticleTable::clusterZMin[A];
123 } 63 }
124 std::fill(candidateConfiguration, << 64 }
125 candidateConfiguration + Parti << 65
126 static_cast<Particle*>(NULL)); << 66 void cleanUp() {
127 << 67 delete candidateConfiguration;
128 std::fill(runningEnergies, << 68 consideredPartners.clear();
129 runningEnergies + ParticleTabl << 69 runningConfiguration.clear();
130 0.0); << 70 }
131 <<
132 std::fill(runningPotentials, <<
133 runningPotentials + ParticleTa <<
134 0.0); <<
135 <<
136 std::fill(runningConfiguration, <<
137 runningConfiguration + Particl <<
138 -1); <<
139 71
>> 72 void zeroOut() {
>> 73 candidateConfiguration = 0;
140 } 74 }
141 75
142 virtual ~ClusteringModelIntercomparison() 76 virtual ~ClusteringModelIntercomparison() {
143 delete [] consideredPartners; << 77 cleanUp();
144 delete [] isInRunningConfiguration; <<
145 } 78 }
146 79
147 virtual Cluster* getCluster(Nucleus*, Part 80 virtual Cluster* getCluster(Nucleus*, Particle*);
148 virtual G4bool clusterCanEscape(Nucleus co << 81 virtual G4bool clusterCanEscape(Cluster const * const);
149 82
150 private: 83 private:
151 void findClusterStartingFrom(const G4int o << 84 void findClusterStartingFrom(const G4int oldA, const G4int oldZ);
152 G4double getPhaseSpace(const G4int oldA, C << 85 G4double getPhaseSpace(G4int oldA, Particle *p);
153 86
154 Nucleus *theNucleus; 87 Nucleus *theNucleus;
>> 88 Particle *theLeadingParticle;
>> 89 ParticleList consideredPartners;
>> 90 ParticleList* candidateConfiguration;
>> 91
>> 92 G4double runningEnergies[ParticleTable::maxClusterMass];
>> 93 ThreeVector runningMomenta[ParticleTable::maxClusterMass];
>> 94 ThreeVector runningPositions[ParticleTable::maxClusterMass];
>> 95 ParticleList runningConfiguration; // Use deque instead?
>> 96 G4double runningPotentials[ParticleTable::maxClusterMass];
155 97
156 G4double runningEnergies[ParticleTable::ma << 98 G4int selectedA, selectedZ;
157 ThreeVector runningMomenta[ParticleTable:: <<
158 ThreeVector runningPositions[ParticleTable <<
159 G4double runningPotentials[ParticleTable:: <<
160 #if defined(INCL_CACHING_CLUSTERING_MODEL_INTE <<
161 Hashing::NucleonItem runningConfiguration[ <<
162 #elif defined(INCL_CACHING_CLUSTERING_MODEL_IN <<
163 G4int runningConfiguration[ParticleTable:: <<
164 #else <<
165 #error Unrecognized INCL_CACHING_CLUSTERING_MO <<
166 #endif <<
167 <<
168 G4int selectedA, selectedZ, selectedS; <<
169 G4double sqtot; 99 G4double sqtot;
170 100
171 G4int clusterZMaxAll, clusterNMaxAll; 101 G4int clusterZMaxAll, clusterNMaxAll;
172 102
173 G4double cascadingEnergyPool; << 103 G4double participantEnergyPool;
174 <<
175 /// \brief Lower limit of Z for cluster of <<
176 static const G4int clusterZMin[ParticleTab <<
177 /// \brief Upper limit of Z for cluster of <<
178 static const G4int clusterZMax[ParticleTab <<
179 <<
180 /// \brief Precomputed factor 1.0/A <<
181 static const G4double clusterPosFact[Parti <<
182 <<
183 /// \brief Precomputed factor (1.0/A)^2 <<
184 static const G4double clusterPosFact2[Part <<
185 <<
186 /// \brief Phase-space parameters for clus <<
187 static const G4double clusterPhaseSpaceCut <<
188 104
189 static const G4double limitCosEscapeAngle; 105 static const G4double limitCosEscapeAngle;
190 <<
191 const G4double protonMass; <<
192 const G4double neutronMass; <<
193 const G4double lambdaMass; <<
194 <<
195 G4int runningMaxClusterAlgorithmMass; <<
196 <<
197 G4int nConsideredMax; <<
198 G4int nConsidered; <<
199 <<
200 /** \brief Array of considered cluster par <<
201 * <<
202 * A dynamical array of ConsideredPartner <<
203 * variable and filled with pointers to nu <<
204 * clustering. We used to use a ParticleLi <<
205 * made it very cumbersome to check whethe <<
206 * included in the running configuration. <<
207 * coupled with a boolean mask (\see{isInR <<
208 * overhead by a large amount. Running ti <<
209 * by almost 30% (!). <<
210 * <<
211 * Lesson learnt: when you need speed, not <<
212 */ <<
213 ConsideredPartner *consideredPartners; <<
214 <<
215 /** \brief Array of flags for nucleons in <<
216 * <<
217 * Clustering partners that are already us <<
218 * configuration are flagged as "true" in <<
219 */ <<
220 G4bool *isInRunningConfiguration; <<
221 <<
222 /** \brief Best cluster configuration <<
223 * <<
224 * This array contains pointers to the nuc <<
225 * cluster configuration that has been fou <<
226 */ <<
227 Particle *candidateConfiguration[ParticleT <<
228 <<
229 #if defined(INCL_CACHING_CLUSTERING_MODEL_INTE <<
230 typedef std::set<Hashing::HashType> HashCo <<
231 typedef HashContainer::iterator HashIterat <<
232 <<
233 /// \brief Array of containers for configu <<
234 HashContainer checkedConfigurations[Partic <<
235 #elif defined(INCL_CACHING_CLUSTERING_MODEL_IN <<
236 /** \brief Class for storing and comparing <<
237 * <<
238 * This class is actually just a wrapper a <<
239 * pointers. It provides a lexicographical <<
240 * (SortedNucleonConfiguration::operator<) <<
241 * containers. <<
242 */ <<
243 class SortedNucleonConfiguration { <<
244 public: <<
245 // Use Particle* as nucleon identifier <<
246 typedef G4int NucleonItem; <<
247 <<
248 /// \brief Constructor <<
249 SortedNucleonConfiguration() : theSize <<
250 <<
251 /// \brief Copy constructor <<
252 SortedNucleonConfiguration(const Sorte <<
253 theSize(rhs.theSize), <<
254 nucleons(new NucleonItem[theSize]) <<
255 { <<
256 std::copy(rhs.nucleons, rhs.nucleons+t <<
257 } <<
258 <<
259 /// \brief Destructor <<
260 ~SortedNucleonConfiguration() { <<
261 delete [] nucleons; <<
262 } <<
263 <<
264 /// \brief Helper method for the assig <<
265 void swap(SortedNucleonConfiguration & <<
266 std::swap(theSize, rhs.theSize); <<
267 std::swap(nucleons, rhs.nucleons); <<
268 } <<
269 <<
270 /// \brief Assignment operator <<
271 SortedNucleonConfiguration &operator=( <<
272 SortedNucleonConfiguration tempConfi <<
273 swap(tempConfig); <<
274 return *this; <<
275 } <<
276 <<
277 /** \brief Order operator for SortedNu <<
278 * <<
279 * The comparison is done lexicographi <<
280 * element to the last). <<
281 */ <<
282 G4bool operator<(const SortedNucleonCo <<
283 // assert(theSize==rhs.theSize); <<
284 return std::lexicographical_compare( <<
285 } <<
286 <<
287 /// \brief Fill configuration with arr <<
288 void fill(NucleonItem *config, size_t <<
289 theSize = n; <<
290 nucleons = new NucleonItem[theSize]; <<
291 std::copy(config, config+theSize, nu <<
292 std::sort(nucleons, nucleons+theSize <<
293 } <<
294 <<
295 private: <<
296 /// \brief Size of the array <<
297 size_t theSize; <<
298 <<
299 /// \brief The real array <<
300 NucleonItem *nucleons; <<
301 }; <<
302 <<
303 typedef std::set<SortedNucleonConfiguratio <<
304 typedef SortedNucleonConfigurationContaine <<
305 <<
306 /// \brief Array of containers for configu <<
307 SortedNucleonConfigurationContainer checke <<
308 #elif !defined(INCL_CACHING_CLUSTERING_MODEL_I <<
309 #error Unrecognized INCL_CACHING_CLUSTERING_MO <<
310 #endif <<
311 <<
312 /** \brief Maximum mass for configuration <<
313 * <<
314 * Skipping configurations becomes ineffic <<
315 */ <<
316 G4int maxMassConfigurationSkipping; <<
317 }; 106 };
318 107
319 } 108 }
320 109
321 #endif 110 #endif
322 111