Geant4 Cross Reference

Cross-Referencing   Geant4
Geant4/processes/hadronic/models/inclxx/incl_physics/include/G4INCLClusteringModelIntercomparison.hh

Version: [ ReleaseNotes ] [ 1.0 ] [ 1.1 ] [ 2.0 ] [ 3.0 ] [ 3.1 ] [ 3.2 ] [ 4.0 ] [ 4.0.p1 ] [ 4.0.p2 ] [ 4.1 ] [ 4.1.p1 ] [ 5.0 ] [ 5.0.p1 ] [ 5.1 ] [ 5.1.p1 ] [ 5.2 ] [ 5.2.p1 ] [ 5.2.p2 ] [ 6.0 ] [ 6.0.p1 ] [ 6.1 ] [ 6.2 ] [ 6.2.p1 ] [ 6.2.p2 ] [ 7.0 ] [ 7.0.p1 ] [ 7.1 ] [ 7.1.p1 ] [ 8.0 ] [ 8.0.p1 ] [ 8.1 ] [ 8.1.p1 ] [ 8.1.p2 ] [ 8.2 ] [ 8.2.p1 ] [ 8.3 ] [ 8.3.p1 ] [ 8.3.p2 ] [ 9.0 ] [ 9.0.p1 ] [ 9.0.p2 ] [ 9.1 ] [ 9.1.p1 ] [ 9.1.p2 ] [ 9.1.p3 ] [ 9.2 ] [ 9.2.p1 ] [ 9.2.p2 ] [ 9.2.p3 ] [ 9.2.p4 ] [ 9.3 ] [ 9.3.p1 ] [ 9.3.p2 ] [ 9.4 ] [ 9.4.p1 ] [ 9.4.p2 ] [ 9.4.p3 ] [ 9.4.p4 ] [ 9.5 ] [ 9.5.p1 ] [ 9.5.p2 ] [ 9.6 ] [ 9.6.p1 ] [ 9.6.p2 ] [ 9.6.p3 ] [ 9.6.p4 ] [ 10.0 ] [ 10.0.p1 ] [ 10.0.p2 ] [ 10.0.p3 ] [ 10.0.p4 ] [ 10.1 ] [ 10.1.p1 ] [ 10.1.p2 ] [ 10.1.p3 ] [ 10.2 ] [ 10.2.p1 ] [ 10.2.p2 ] [ 10.2.p3 ] [ 10.3 ] [ 10.3.p1 ] [ 10.3.p2 ] [ 10.3.p3 ] [ 10.4 ] [ 10.4.p1 ] [ 10.4.p2 ] [ 10.4.p3 ] [ 10.5 ] [ 10.5.p1 ] [ 10.6 ] [ 10.6.p1 ] [ 10.6.p2 ] [ 10.6.p3 ] [ 10.7 ] [ 10.7.p1 ] [ 10.7.p2 ] [ 10.7.p3 ] [ 10.7.p4 ] [ 11.0 ] [ 11.0.p1 ] [ 11.0.p2 ] [ 11.0.p3, ] [ 11.0.p4 ] [ 11.1 ] [ 11.1.1 ] [ 11.1.2 ] [ 11.1.3 ] [ 11.2 ] [ 11.2.1 ] [ 11.2.2 ] [ 11.3.0 ]

  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 // INCL++ intra-nuclear cascade model
 27 // Alain Boudard, CEA-Saclay, France
 28 // Joseph Cugnon, University of Liege, Belgium
 29 // Jean-Christophe David, CEA-Saclay, France
 30 // Pekka Kaitaniemi, CEA-Saclay, France, and Helsinki Institute of Physics, Finland
 31 // Sylvie Leray, CEA-Saclay, France
 32 // Davide Mancusi, CEA-Saclay, France
 33 //
 34 #define INCLXX_IN_GEANT4_MODE 1
 35 
 36 #include "globals.hh"
 37 
 38 #ifndef G4INCLClusteringModelIntercomparison_hh
 39 #define G4INCLClusteringModelIntercomparison_hh 1
 40 
 41 #ifdef INCLXX_IN_GEANT4_MODE
 42 #define INCL_CACHING_CLUSTERING_MODEL_INTERCOMPARISON_Set 1
 43 #endif // INCLXX_IN_GEANT4_MODE
 44 
 45 #include "G4INCLIClusteringModel.hh"
 46 #include "G4INCLParticle.hh"
 47 #include "G4INCLParticleTable.hh"
 48 #include "G4INCLCluster.hh"
 49 #include "G4INCLNucleus.hh"
 50 #include "G4INCLKinematicsUtils.hh"
 51 #include "G4INCLHashing.hh"
 52 
 53 #include <set>
 54 #include <algorithm>
 55 
 56 namespace G4INCL {
 57 
 58   /** \brief Container for the relevant information
 59    *
 60    * This struct contains all the information that is relevant for the
 61    * clustering algorithm. It is probably more compact than the Particles it
 62    * feeds on, hopefully improving cache performance.
 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->isTargetSpectator()),
 86       Z(particle->getZ()),
 87       S(particle->getS()),
 88       position(particle->getPosition()),
 89       momentum(particle->getMomentum()),
 90       energy(particle->getEnergy()),
 91       potentialEnergy(particle->getPotentialEnergy())
 92     {}
 93   };
 94 
 95   /// \brief Cluster coalescence algorithm used in the IAEA intercomparison
 96   class ClusteringModelIntercomparison : public IClusteringModel {
 97   public:
 98     ClusteringModelIntercomparison(Config const * const theConfig) :
 99       theNucleus(NULL),
100       selectedA(0),
101       selectedZ(0),
102       selectedS(0),
103       sqtot(0.),
104       cascadingEnergyPool(0.),
105       protonMass(ParticleTable::getRealMass(Proton)),
106       neutronMass(ParticleTable::getRealMass(Neutron)),
107       lambdaMass(ParticleTable::getRealMass(Lambda)),
108       runningMaxClusterAlgorithmMass(theConfig->getClusterMaxMass()),
109       nConsideredMax(0),
110       nConsidered(0),
111       consideredPartners(NULL),
112       isInRunningConfiguration(NULL),
113       maxMassConfigurationSkipping(ParticleTable::maxClusterMass)
114     {
115       // Set up the maximum charge and neutron number for clusters
116       clusterZMaxAll = 0;
117       clusterNMaxAll = 0;
118       for(G4int A=0; A<=runningMaxClusterAlgorithmMass; ++A) {
119         if(clusterZMax[A]>clusterZMaxAll)
120           clusterZMaxAll = clusterZMax[A];
121         if(A-clusterZMin[A]>clusterNMaxAll)
122           clusterNMaxAll = A-clusterZMin[A];
123       }
124       std::fill(candidateConfiguration,
125                 candidateConfiguration + ParticleTable::maxClusterMass,
126                 static_cast<Particle*>(NULL));
127 
128       std::fill(runningEnergies,
129                 runningEnergies + ParticleTable::maxClusterMass,
130                 0.0);
131 
132       std::fill(runningPotentials,
133                 runningPotentials + ParticleTable::maxClusterMass,
134                 0.0);
135 
136       std::fill(runningConfiguration,
137                 runningConfiguration + ParticleTable::maxClusterMass,
138                 -1);
139 
140     }
141 
142     virtual ~ClusteringModelIntercomparison() {
143       delete [] consideredPartners;
144       delete [] isInRunningConfiguration;
145     }
146 
147     virtual Cluster* getCluster(Nucleus*, Particle*);
148     virtual G4bool clusterCanEscape(Nucleus const * const, Cluster const * const);
149 
150   private:
151     void findClusterStartingFrom(const G4int oldA, const G4int oldZ, const G4int oldS);
152     G4double getPhaseSpace(const G4int oldA, ConsideredPartner const &p);
153 
154     Nucleus *theNucleus;
155 
156     G4double runningEnergies[ParticleTable::maxClusterMass+1];
157     ThreeVector runningMomenta[ParticleTable::maxClusterMass+1];
158     ThreeVector runningPositions[ParticleTable::maxClusterMass+1];
159     G4double runningPotentials[ParticleTable::maxClusterMass+1];
160 #if defined(INCL_CACHING_CLUSTERING_MODEL_INTERCOMPARISON_HashMask)
161     Hashing::NucleonItem runningConfiguration[ParticleTable::maxClusterMass];
162 #elif defined(INCL_CACHING_CLUSTERING_MODEL_INTERCOMPARISON_Set) || defined(INCL_CACHING_CLUSTERING_MODEL_INTERCOMPARISON_None)
163     G4int runningConfiguration[ParticleTable::maxClusterMass];
164 #else
165 #error Unrecognized INCL_CACHING_CLUSTERING_MODEL_INTERCOMPARISON. Allowed values are: Set, HashMask, None.
166 #endif
167 
168     G4int selectedA, selectedZ, selectedS;
169     G4double sqtot;
170 
171     G4int clusterZMaxAll, clusterNMaxAll;
172 
173     G4double cascadingEnergyPool;
174 
175     /// \brief Lower limit of Z for cluster of mass A
176     static const G4int clusterZMin[ParticleTable::maxClusterMass+1];
177     /// \brief Upper limit of Z for cluster of mass A
178     static const G4int clusterZMax[ParticleTable::maxClusterMass+1];
179 
180     /// \brief Precomputed factor 1.0/A
181     static const G4double clusterPosFact[ParticleTable::maxClusterMass+1];
182 
183     /// \brief Precomputed factor (1.0/A)^2
184     static const G4double clusterPosFact2[ParticleTable::maxClusterMass+1];
185 
186     /// \brief Phase-space parameters for cluster formation
187     static const G4double clusterPhaseSpaceCut[ParticleTable::maxClusterMass+1];
188 
189     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 partners
201      *
202      * A dynamical array of ConsideredPartner objects is allocated on this
203      * variable and filled with pointers to nucleons which are eligible for
204      * clustering. We used to use a ParticleList for this purpose, but this
205      * made it very cumbersome to check whether nucleons had already been
206      * included in the running configuration. Using an array of Particle*
207      * coupled with a boolean mask (\see{isInRunningConfiguration}) reduces the
208      * overhead by a large amount.  Running times for 1-GeV p+Pb208 went down
209      * by almost 30% (!).
210      *
211      * Lesson learnt: when you need speed, nothing beats a good ol' array.
212      */
213     ConsideredPartner *consideredPartners;
214 
215     /** \brief Array of flags for nucleons in the running configuration
216      *
217      * Clustering partners that are already used in the running cluster
218      * configuration are flagged as "true" in this array.
219      */
220     G4bool *isInRunningConfiguration;
221 
222     /** \brief Best cluster configuration
223      *
224      * This array contains pointers to the nucleons which make up the best
225      * cluster configuration that has been found so far.
226      */
227     Particle *candidateConfiguration[ParticleTable::maxClusterMass];
228 
229 #if defined(INCL_CACHING_CLUSTERING_MODEL_INTERCOMPARISON_HashMask)
230     typedef std::set<Hashing::HashType> HashContainer;
231     typedef HashContainer::iterator HashIterator;
232 
233     /// \brief Array of containers for configurations that have already been checked
234     HashContainer checkedConfigurations[ParticleTable::maxClusterMass-2];
235 #elif defined(INCL_CACHING_CLUSTERING_MODEL_INTERCOMPARISON_Set)
236     /** \brief Class for storing and comparing sorted nucleon configurations
237      *
238      * This class is actually just a wrapper around an array of Particle*
239      * pointers. It provides a lexicographical comparison operator
240      * (SortedNucleonConfiguration::operator<) for inclusion in std::set
241      * containers.
242      */
243     class SortedNucleonConfiguration {
244       public:
245         // Use Particle* as nucleon identifiers
246         typedef G4int NucleonItem;
247 
248         /// \brief Constructor
249         SortedNucleonConfiguration() : theSize(0), nucleons(NULL) {}
250 
251         /// \brief Copy constructor
252         SortedNucleonConfiguration(const SortedNucleonConfiguration &rhs) :
253           theSize(rhs.theSize),
254           nucleons(new NucleonItem[theSize])
255       {
256         std::copy(rhs.nucleons, rhs.nucleons+theSize, nucleons);
257       }
258 
259         /// \brief Destructor
260         ~SortedNucleonConfiguration() {
261           delete [] nucleons;
262         }
263 
264         /// \brief Helper method for the assignment operator
265         void swap(SortedNucleonConfiguration &rhs) {
266           std::swap(theSize, rhs.theSize);
267           std::swap(nucleons, rhs.nucleons);
268         }
269 
270         /// \brief Assignment operator
271         SortedNucleonConfiguration &operator=(const SortedNucleonConfiguration &rhs) {
272           SortedNucleonConfiguration tempConfig(rhs);
273           swap(tempConfig);
274           return *this;
275         }
276 
277         /** \brief Order operator for SortedNucleonConfiguration
278          *
279          * The comparison is done lexicographically (i.e. from the first
280          * element to the last).
281          */
282         G4bool operator<(const SortedNucleonConfiguration &rhs) const {
283 // assert(theSize==rhs.theSize);
284           return std::lexicographical_compare(nucleons, nucleons+theSize, rhs.nucleons, rhs.nucleons+theSize);
285         }
286 
287         /// \brief Fill configuration with array of NucleonItem
288         void fill(NucleonItem *config, size_t n) {
289           theSize = n;
290           nucleons = new NucleonItem[theSize];
291           std::copy(config, config+theSize, nucleons);
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<SortedNucleonConfiguration> SortedNucleonConfigurationContainer;
304     typedef SortedNucleonConfigurationContainer::iterator SortedNucleonConfigurationIterator;
305 
306     /// \brief Array of containers for configurations that have already been checked
307     SortedNucleonConfigurationContainer checkedConfigurations[ParticleTable::maxClusterMass-2];
308 #elif !defined(INCL_CACHING_CLUSTERING_MODEL_INTERCOMPARISON_None)
309 #error Unrecognized INCL_CACHING_CLUSTERING_MODEL_INTERCOMPARISON. Allowed values are: Set, HashMask, None.
310 #endif
311 
312     /** \brief Maximum mass for configuration storage
313      *
314      * Skipping configurations becomes inefficient above this mass.
315      */
316     G4int maxMassConfigurationSkipping;
317   };
318 
319 }
320 
321 #endif
322