Geant4  10.02.p01
G4Evaporation.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 //
27 // $Id: G4Evaporation.cc 94459 2015-11-18 14:39:46Z gcosmo $
28 //
29 // Hadronic Process: Nuclear De-excitations
30 // by V. Lara (Oct 1998)
31 //
32 // Alex Howard - added protection for negative probabilities in the sum, 14/2/07
33 //
34 // Modif (03 September 2008) by J. M. Quesada for external choice of inverse
35 // cross section option
36 // JMQ (06 September 2008) Also external choices have been added for
37 // superimposed Coulomb barrier (if useSICBis set true, by default is false)
38 //
39 // V.Ivanchenko (27 July 2009) added G4EvaporationDefaultGEMFactory option
40 // V.Ivanchenko (10 May 2010) rewrited BreakItUp method: do not make new/delete
41 // photon channel first, fission second,
42 // added G4UnstableFragmentBreakUp to decay
43 // unphysical fragments (like 2n or 2p),
44 // use Z and A integer
45 // V.Ivanchenko (22 April 2011) added check if a fragment can be deexcited
46 // by the FermiBreakUp model
47 // V.Ivanchenko (23 January 2012) added pointer of G4VPhotonEvaporation
48 // V.Ivanchenko (6 May 2013) added check of existence of residual ion
49 // in the ion table
50 
51 #include "G4Evaporation.hh"
52 #include "G4SystemOfUnits.hh"
53 #include "G4EvaporationFactory.hh"
56 #include "G4HadronicException.hh"
57 #include "G4NistManager.hh"
58 #include "G4FermiFragmentsPool.hh"
59 #include "G4PhotonEvaporation.hh"
61 #include "G4VEvaporationChannel.hh"
62 #include "G4ParticleTable.hh"
63 #include "G4IonTable.hh"
64 
66  : nChannels(0)
67 {
69 
70  G4VEvaporationChannel* gammaEvap = 0;
71  char* env = getenv("G4UsePhotonEvaporationOLD");
72  if(env) { gammaEvap = new G4PhotonEvaporationOLD(); }
73  else { gammaEvap = new G4PhotonEvaporation(); }
74  SetPhotonEvaporation(gammaEvap);
75 
77  SetParameters();
80 }
81 
83  : nChannels(0)
84 {
85  if(photoEvaporation) { SetPhotonEvaporation(photoEvaporation); }
87 
90  SetParameters();
93 }
94 
96 {}
97 
99 {
104  probabilities.reserve(68);
105 }
106 
108 {
109  CleanChannels();
111  nChannels = theChannels->size();
112  probabilities.resize(nChannels, 0.0);
114 }
115 
117 {
118  for(size_t i=0; i<nChannels; ++i) {
119  (*theChannels)[i]->SetOPTxs(OPTxs);
120  (*theChannels)[i]->UseSICB(useSICB);
121  (*theChannels)[i]->Initialise();
122  }
123 }
124 
126 {
127  delete theChannelFactory;
130 }
131 
133 {
134  delete theChannelFactory;
137 }
138 
140 {
141  delete theChannelFactory;
144 }
145 
147 {
148  G4FragmentVector * theResult = new G4FragmentVector;
149  G4Fragment* theResidualNucleus = new G4Fragment(theNucleus);
150  BreakFragment(theResult, theResidualNucleus);
151  return theResult;
152 }
153 
155  G4Fragment* theResidualNucleus)
156 {
157  G4double totprob, prob, oldprob = 0.0;
158  size_t maxchannel, i;
159 
160  G4int Amax = theResidualNucleus->GetA_asInt();
161 
162  // Starts loop over evaporated particles, loop is limited by number
163  // of nucleons
164  for(G4int ia=0; ia<Amax; ++ia) {
165 
166  // g,n,p and light fragments - evaporation is finished
167  G4int Z = theResidualNucleus->GetZ_asInt();
168  G4int A = theResidualNucleus->GetA_asInt();
169  G4double Eex = theResidualNucleus->GetExcitationEnergy();
170  G4double mass = theResidualNucleus->GetGroundStateMass();
171 
172  // stop deecitation loop if residual can be deexcited by FBU
173  if(maxZforFBU > Z && maxAforFBU > A && Z > 0 && A > Z) {
174  if(thePool->IsApplicable(Z, A, mass+Eex)) {
175  theResult->push_back(theResidualNucleus);
176  return;
177  }
178  }
179  // check if it is stable, then finish evaporation
180  G4double abun = nist->GetIsotopeAbundance(Z, A);
181  /*
182  G4cout << "### G4Evaporation::BreakItUp step " << ia << " Z= " << Z
183  << " A= " << A << " Eex(MeV)= "
184  << theResidualNucleus->GetExcitationEnergy()
185  << " aban= " << abun << G4endl;
186  */
187  // stop deecitation loop in the case of a cold stable fragment
188  if(Eex <= minExcitation && abun > 0.0) {
189  theResult->push_back(theResidualNucleus);
190  return;
191  }
192 
193  totprob = 0.0;
194  maxchannel = nChannels;
195  /*
196  G4cout << "### Evaporation loop #" << ia
197  << " Fragment: " << theResidualNucleus << G4endl;
198  */
199  // loop over evaporation channels
200  for(i=0; i<nChannels; ++i) {
201  prob = (*theChannels)[i]->GetEmissionProbability(theResidualNucleus);
202  // G4cout << " Channel# " << i << " prob= " << prob << G4endl;
203 
204  totprob += prob;
205  probabilities[i] = totprob;
206 
207  // if two recent probabilities are near zero stop computations
208  if(i>=8) {
209  if(prob <= totprob*1.e-8 && oldprob <= totprob*1.e-8) {
210  maxchannel = i+1;
211  break;
212  }
213  }
214  oldprob = prob;
215  }
216 
217  // photon evaporation in the case of no other channels available
218  // do evaporation chain and reset total probability
219  if(0.0 < totprob && probabilities[0] == totprob) {
220  // G4cout << "Start chain of gamma evaporation" << G4endl;
221  (*theChannels)[0]->BreakUpChain(theResult, theResidualNucleus);
222  totprob = 0.0;
223  }
224 
225  // stable fragment - evaporation is finished
226  if(0.0 == totprob) {
227 
228  // if fragment is exotic, then force its decay
229  if(0.0 == abun && Z < 20) {
230  //G4cout << "$$$ Decay exotic fragment" << G4endl;
231  unstableBreakUp.BreakUpChain(theResult, theResidualNucleus);
232  } else {
233  theResult->push_back(theResidualNucleus);
234  }
235  return;
236  }
237 
238  // select channel
239  totprob *= G4UniformRand();
240  // loop over evaporation channels
241  for(i=0; i<maxchannel; ++i) { if(probabilities[i] >= totprob) { break; } }
242 
243  //G4cout << "Channel # " << i << G4endl;
244  G4Fragment* frag = (*theChannels)[i]->EmittedFragment(theResidualNucleus);
245  //if(frag) G4cout << " " << *frag << G4endl;
246  if(frag) { theResult->push_back(frag); }
247  // selected channel cannot sample secondary
248  else {
249  theResult->push_back(theResidualNucleus);
250  return;
251  }
252  }
253 
254  // loop is stopped, save residual, which is unclear state
255  theResult->push_back(theResidualNucleus);
256 }
G4FragmentVector * BreakItUp(const G4Fragment &theNucleus)
void BreakFragment(G4FragmentVector *, G4Fragment *theNucleus)
void SetDefaultChannel()
virtual std::vector< G4VEvaporationChannel * > * GetChannel()=0
std::vector< G4double > probabilities
G4VEvaporationChannel * thePhotonEvaporation
int G4int
Definition: G4Types.hh:78
static G4NistManager * Instance()
std::vector< G4VEvaporationChannel * > * theChannels
void SetGEMChannel()
virtual void InitialiseChannels()
G4bool IsApplicable(G4int Z, G4int A, G4double mass) const
G4IonTable * GetIonTable() const
#define G4UniformRand()
Definition: Randomize.hh:97
double A(double temperature)
G4int GetA_asInt() const
Definition: G4Fragment.hh:251
std::vector< G4Fragment * > G4FragmentVector
Definition: G4Fragment.hh:63
G4IonTable * theTableOfIons
G4double GetIsotopeAbundance(G4int Z, G4int N) const
G4double GetGroundStateMass() const
Definition: G4Fragment.hh:273
G4NistManager * nist
void SetCombinedChannel()
G4FermiFragmentsPool * thePool
static G4ParticleTable * GetParticleTable()
G4double minExcitation
void InitialiseEvaporation()
virtual void SetPhotonEvaporation(G4VEvaporationChannel *ptr)
G4UnstableFragmentBreakUp unstableBreakUp
virtual ~G4Evaporation()
G4int GetZ_asInt() const
Definition: G4Fragment.hh:256
void SetParameters()
static const double keV
Definition: G4SIunits.hh:213
double G4double
Definition: G4Types.hh:76
virtual G4bool BreakUpChain(G4FragmentVector *theResult, G4Fragment *theNucleus)
static G4FermiFragmentsPool * Instance()
G4VEvaporationFactory * theChannelFactory
G4double GetExcitationEnergy() const
Definition: G4Fragment.hh:268