Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
•All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros
G4FermiFragmentsPool.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// $Id: G4VFermiBreakUp.cc,v 1.5 2006-06-29 20:13:13 gunter Exp $
27// GEANT4 tag $Name: not supported by cvs2svn $
28//
29// Hadronic Process: Nuclear De-excitations
30// by V. Lara
31//
32// Modifications:
33// J.M.Quesada, July 2009, bug fixed in excitation energies:
34// ALL of them are in MeV instead of keV (as they were expressed previously)
35// source: http://www.nndc.bnl.gov/chart
36// Unknown excitation energies in He5 and Li5 have been suppressed
37// Long lived levels (half-lives of the order ps-fs have been included)
38//
39// J. M. Quesada, April 2010: excitation energies according to tabulated values
40// in PhotonEvaporatoion2.0. Fake photons eliminated.
41//
42// 01.04.2011 General cleanup by V.Ivanchenko - more clean usage of static
43//
44// 04.05.2011 J. M. Quesada: added detailed printout for testing
45
46#include "G4FermiFragmentsPool.hh"
47#include "G4PhysicalConstants.hh"
48#include "G4SystemOfUnits.hh"
49#include "G4StableFermiFragment.hh"
50#include "G4B9FermiFragment.hh"
51#include "G4Be8FermiFragment.hh"
52#include "G4He5FermiFragment.hh"
53#include "G4Li5FermiFragment.hh"
54
55G4FermiFragmentsPool* G4FermiFragmentsPool::theInstance = 0;
56
57G4FermiFragmentsPool* G4FermiFragmentsPool::Instance()
58{
59 if(0 == theInstance) {
60 static G4FermiFragmentsPool pool;
61 theInstance = &pool;
62 }
63 return theInstance;
64}
65
66G4FermiFragmentsPool::G4FermiFragmentsPool()
67{
68 maxZ = 9;
69 maxA = 17;
70 verbose = 0;
71 Initialise();
72}
73
74G4FermiFragmentsPool::~G4FermiFragmentsPool()
75{
76 for(size_t i=0; i<17; ++i) {
77 size_t nn = list1[i].size();
78 if(0 < nn) { for(size_t j=0; j<nn; ++j) { delete (list1[i])[j]; }}
79 nn = list2[i].size();
80 if(0 < nn) { for(size_t j=0; j<nn; ++j) { delete (list2[i])[j]; }}
81 nn = list3[i].size();
82 if(0 < nn) { for(size_t j=0; j<nn; ++j) { delete (list3[i])[j]; }}
83 nn = list4[i].size();
84 if(0 < nn) { for(size_t j=0; j<nn; ++j) { delete (list4[i])[j]; }}
85 }
86 size_t nn = listextra.size();
87 if(0 < nn) { for(size_t j=0; j<nn; ++j) { delete listextra[j]; }}
88 nn = fragment_pool.size();
89 if(0 < nn) { for(size_t j=0; j<nn; ++j) { delete fragment_pool[j]; }}
90}
91
92void G4FermiFragmentsPool::Initialise()
93{
94 // JMQ 30/06/09 unknown levels have been supressed
95 // JMQ 01/07/09 corrected excitation energies for 64-66, according to
96 // http://www.nndc.bnl.gov/chart
97 // JMQ 19/04/10 new level, fragment numbering shifted accordingly from here onwards
98 // A Z Pol ExcitE
99 fragment_pool.push_back(new G4StableFermiFragment( 1, 0, 2, 0.00*MeV ));
100 fragment_pool.push_back(new G4StableFermiFragment( 1, 1, 2, 0.00*MeV ));
101 fragment_pool.push_back(new G4StableFermiFragment( 2, 1, 3, 0.00*MeV ));
102 fragment_pool.push_back(new G4StableFermiFragment( 3, 1, 2, 0.00*MeV ));
103 fragment_pool.push_back(new G4StableFermiFragment( 3, 2, 2, 0.00*MeV ));
104 fragment_pool.push_back(new G4StableFermiFragment( 4, 2, 1, 0.00*MeV ));
105 fragment_pool.push_back(new G4He5FermiFragment ( 5, 2, 4, 0.00*MeV ));
106 fragment_pool.push_back(new G4Li5FermiFragment ( 5, 3, 4, 0.00*MeV ));
107 fragment_pool.push_back(new G4StableFermiFragment( 6, 2, 1, 0.00*MeV ));
108 fragment_pool.push_back(new G4StableFermiFragment( 6, 3, 3, 0.00*MeV ));
109 fragment_pool.push_back(new G4StableFermiFragment( 6, 3, 1, 3.562880*MeV ));
110 fragment_pool.push_back(new G4StableFermiFragment( 7, 3, 4, 0.00*MeV ));
111 fragment_pool.push_back(new G4StableFermiFragment( 7, 3, 2, 0.4776120*MeV ));
112 fragment_pool.push_back(new G4StableFermiFragment( 7, 4, 4, 0.00*MeV ));
113 fragment_pool.push_back(new G4StableFermiFragment( 7, 4, 2, 0.4290800*MeV ));
114 fragment_pool.push_back(new G4StableFermiFragment( 8, 3, 5, 0.00*MeV ));
115 fragment_pool.push_back(new G4StableFermiFragment( 8, 3, 3, 0.9808000*MeV ));
116 fragment_pool.push_back(new G4Be8FermiFragment ( 8, 4, 1, 0.00*MeV ));
117 fragment_pool.push_back(new G4StableFermiFragment( 9, 4, 4, 0.00*MeV ));
118 fragment_pool.push_back(new G4B9FermiFragment ( 9, 5, 4, 0.00*MeV ));
119 fragment_pool.push_back(new G4StableFermiFragment( 10, 4, 1, 0.00*MeV ));
120 fragment_pool.push_back(new G4StableFermiFragment( 10, 4, 5, 3.368030*MeV ));
121 fragment_pool.push_back(new G4StableFermiFragment( 10, 4, 8, 5.958390*MeV ));
122 fragment_pool.push_back(new G4StableFermiFragment( 10, 4, 1, 6.179300*MeV ));
123 fragment_pool.push_back(new G4StableFermiFragment( 10, 4, 5, 6.263300*MeV ));
124 fragment_pool.push_back(new G4StableFermiFragment( 10, 5, 7, 0.00*MeV ));
125 fragment_pool.push_back(new G4StableFermiFragment( 10, 5, 3, 0.7183500*MeV ));
126 fragment_pool.push_back(new G4StableFermiFragment( 10, 5, 1, 1.740150*MeV ));
127 fragment_pool.push_back(new G4StableFermiFragment( 10, 5, 3, 2.154300*MeV ));
128 fragment_pool.push_back(new G4StableFermiFragment( 10, 5, 5, 3.587100*MeV ));
129 fragment_pool.push_back(new G4StableFermiFragment( 10, 6, 3, 0.00*MeV ));
130 fragment_pool.push_back(new G4StableFermiFragment( 10, 6, 5, 3.353600*MeV ));
131 fragment_pool.push_back(new G4StableFermiFragment( 11, 5, 4, 0.00*MeV ));
132 fragment_pool.push_back(new G4StableFermiFragment( 11, 5, 2, 2.124693*MeV ));
133 fragment_pool.push_back(new G4StableFermiFragment( 11, 5, 6, 4.444890*MeV ));
134 fragment_pool.push_back(new G4StableFermiFragment( 11, 5, 4, 5.020310*MeV ));
135 fragment_pool.push_back(new G4StableFermiFragment( 11, 5, 8, 6.742900*MeV ));
136 fragment_pool.push_back(new G4StableFermiFragment( 11, 5, 2, 6.791800*MeV ));
137 fragment_pool.push_back(new G4StableFermiFragment( 11, 5, 6, 7.285510*MeV ));
138 fragment_pool.push_back(new G4StableFermiFragment( 11, 5, 4, 7.977840*MeV ));
139 fragment_pool.push_back(new G4StableFermiFragment( 11, 5, 6, 8.560300*MeV ));
140 fragment_pool.push_back(new G4StableFermiFragment( 11, 6, 4, 0.00*MeV ));
141 fragment_pool.push_back(new G4StableFermiFragment( 11, 6, 2, 2.00*MeV ));
142 fragment_pool.push_back(new G4StableFermiFragment( 11, 6, 6, 4.318800*MeV ));
143 fragment_pool.push_back(new G4StableFermiFragment( 11, 6, 4, 4.804200*MeV ));
144 fragment_pool.push_back(new G4StableFermiFragment( 11, 6, 2, 6.339200*MeV ));
145 fragment_pool.push_back(new G4StableFermiFragment( 11, 6, 8, 6.478200*MeV ));
146 fragment_pool.push_back(new G4StableFermiFragment( 11, 6, 6, 6.904800*MeV ));
147 fragment_pool.push_back(new G4StableFermiFragment( 11, 6, 4, 7.499700*MeV ));
148 fragment_pool.push_back(new G4StableFermiFragment( 11, 6, 4, 8.104500*MeV ));
149 fragment_pool.push_back(new G4StableFermiFragment( 11, 6, 6, 8.420000*MeV ));
150 fragment_pool.push_back(new G4StableFermiFragment( 12, 5, 3, 0.00*MeV ));
151 fragment_pool.push_back(new G4StableFermiFragment( 12, 5, 5, 0.9531400*MeV ));
152 fragment_pool.push_back(new G4StableFermiFragment( 12, 5, 5, 1.673650*MeV ));
153 fragment_pool.push_back(new G4StableFermiFragment( 12, 5, 3, 2.620800*MeV ));
154 fragment_pool.push_back(new G4StableFermiFragment( 12, 6, 1, 0.00*MeV ));
155 fragment_pool.push_back(new G4StableFermiFragment( 12, 6, 5, 4.438910*MeV ));
156 fragment_pool.push_back(new G4StableFermiFragment( 13, 6, 2, 0.00*MeV ));
157 fragment_pool.push_back(new G4StableFermiFragment( 13, 6, 2, 3.089443*MeV ));
158 fragment_pool.push_back(new G4StableFermiFragment( 13, 6, 4, 3.684507*MeV ));
159 fragment_pool.push_back(new G4StableFermiFragment( 13, 6, 6, 3.853807*MeV ));
160 fragment_pool.push_back(new G4StableFermiFragment( 13, 7, 2, 0.00*MeV ));
161 fragment_pool.push_back(new G4StableFermiFragment( 14, 6, 1, 0.00*MeV ));
162 fragment_pool.push_back(new G4StableFermiFragment( 14, 6, 3, 6.093800*MeV ));
163 fragment_pool.push_back(new G4StableFermiFragment( 14, 6, 1, 6.589400*MeV ));
164 fragment_pool.push_back(new G4StableFermiFragment( 14, 6, 7, 6.728200*MeV ));
165 fragment_pool.push_back(new G4StableFermiFragment( 14, 6, 1, 6.902600*MeV ));
166 fragment_pool.push_back(new G4StableFermiFragment( 14, 6, 5, 7.012000*MeV ));
167 fragment_pool.push_back(new G4StableFermiFragment( 14, 6, 5, 7.341000*MeV ));
168 fragment_pool.push_back(new G4StableFermiFragment( 14, 7, 3, 0.00*MeV ));
169 fragment_pool.push_back(new G4StableFermiFragment( 14, 7, 1, 2.312798*MeV ));
170 fragment_pool.push_back(new G4StableFermiFragment( 14, 7, 3, 3.948100*MeV ));
171 fragment_pool.push_back(new G4StableFermiFragment( 14, 7, 1, 4.915100*MeV ));
172 fragment_pool.push_back(new G4StableFermiFragment( 14, 7, 5, 5.105890*MeV ));
173 fragment_pool.push_back(new G4StableFermiFragment( 14, 7, 3, 5.691440*MeV ));
174 fragment_pool.push_back(new G4StableFermiFragment( 14, 7, 7, 5.834250*MeV ));
175 fragment_pool.push_back(new G4StableFermiFragment( 14, 7, 3, 6.203500*MeV ));
176 fragment_pool.push_back(new G4StableFermiFragment( 14, 7, 7, 6.446170*MeV ));
177 fragment_pool.push_back(new G4StableFermiFragment( 14, 7, 5, 7.029120*MeV ));
178 fragment_pool.push_back(new G4StableFermiFragment( 15, 7, 2, 0.00*MeV ));
179 // JMQ 010709 two very close levels instead of only one, with their own spins
180 fragment_pool.push_back(new G4StableFermiFragment( 15, 7, 6, 5.270155*MeV ));
181 fragment_pool.push_back(new G4StableFermiFragment( 15, 7, 2, 5.298822*MeV ));
182 fragment_pool.push_back(new G4StableFermiFragment( 15, 7, 4, 6.323780*MeV ));
183 //JMQ 010709 new level and corrected energy and spins
184 fragment_pool.push_back(new G4StableFermiFragment( 15, 7, 6, 7.155050*MeV ));
185 fragment_pool.push_back(new G4StableFermiFragment( 15, 7, 4, 7.300830*MeV ));
186 fragment_pool.push_back(new G4StableFermiFragment( 15, 7, 8, 7.567100*MeV ));
187 fragment_pool.push_back(new G4StableFermiFragment( 15, 7, 2, 8.312620*MeV ));
188 fragment_pool.push_back(new G4StableFermiFragment( 15, 7, 4, 8.571400*MeV ));
189 fragment_pool.push_back(new G4StableFermiFragment( 15, 7, 2, 9.049710*MeV ));
190 //JMQ 010709 new levels for N15
191 fragment_pool.push_back(new G4StableFermiFragment( 15, 7, 4, 9.151900*MeV ));
192 fragment_pool.push_back(new G4StableFermiFragment( 15, 7, 6, 9.154900*MeV ));
193 fragment_pool.push_back(new G4StableFermiFragment( 15, 7, 2, 9.222100*MeV ));
194 fragment_pool.push_back(new G4StableFermiFragment( 15, 7, 6, 9.760000*MeV ));
195 fragment_pool.push_back(new G4StableFermiFragment( 15, 7, 8, 9.829000*MeV ));
196 fragment_pool.push_back(new G4StableFermiFragment( 15, 7, 4, 9.925000*MeV ));
197 fragment_pool.push_back(new G4StableFermiFragment( 15, 7, 4, 10.06600*MeV ));
198 fragment_pool.push_back(new G4StableFermiFragment( 15, 8, 2, 0.00*MeV ));
199 //JMQ 010709 new level and spins
200 fragment_pool.push_back(new G4StableFermiFragment( 15, 8, 2, 5.183000*MeV ));
201 fragment_pool.push_back(new G4StableFermiFragment( 15, 8, 6, 5.240900*MeV ));
202 fragment_pool.push_back(new G4StableFermiFragment( 15, 8, 4, 6.176300*MeV ));
203 fragment_pool.push_back(new G4StableFermiFragment( 15, 8, 4, 6.793100*MeV ));
204 fragment_pool.push_back(new G4StableFermiFragment( 15, 8, 6, 6.859400*MeV ));
205 fragment_pool.push_back(new G4StableFermiFragment( 15, 8, 8, 7.275900*MeV ));
206 fragment_pool.push_back(new G4StableFermiFragment( 16, 7, 5, 0.00*MeV ));
207 fragment_pool.push_back(new G4StableFermiFragment( 16, 7, 1, 0.1204200*MeV ));
208 fragment_pool.push_back(new G4StableFermiFragment( 16, 7, 7, 0.2982200*MeV ));
209 fragment_pool.push_back(new G4StableFermiFragment( 16, 7, 3, 0.3972700*MeV ));
210 //JMQ 010709 some energies and spins have been changed
211 fragment_pool.push_back(new G4StableFermiFragment( 16, 8, 1, 0.00*MeV ));
212 fragment_pool.push_back(new G4StableFermiFragment( 16, 8, 1, 6.049400*MeV ));
213 fragment_pool.push_back(new G4StableFermiFragment( 16, 8, 7, 6.129890*MeV ));
214 fragment_pool.push_back(new G4StableFermiFragment( 16, 8, 5, 6.917100*MeV ));
215 //JMQ 180510 fixed fragment 111
216 fragment_pool.push_back(new G4StableFermiFragment( 16, 8, 3, 7.116850*MeV ));
217
218 G4int nfrag = fragment_pool.size();
219
220 // list of fragments ordered by A
221 for(G4int i=0; i<nfrag; ++i) {
222 std::vector<const G4VFermiFragment*> newvec;
223 newvec.push_back(fragment_pool[i]);
224 G4FermiConfiguration* conf = new G4FermiConfiguration(newvec);
225 G4int A = fragment_pool[i]->GetA();
226 list1[A].push_back(conf);
227 }
228 if(verbose > 0) {
229 G4cout << "### G4FermiFragmentPool: " << nfrag
230 << " fragments" << G4endl;
231 for(G4int A=1; A<maxA; ++A) {
232 G4cout << " A= " << A << " : Z= ";
233 for(size_t j=0; j<list1[A].size(); ++j) {
234 G4cout << (list1[A])[j]->GetZ() << " ";
235 }
236 G4cout << G4endl;
237 }
238 }
239
240 // list of fragment pairs ordered by A
241 G4int counter = 0;
242 G4int tot = 0;
243 for(G4int i=0; i<nfrag; ++i) {
244 G4int Z1 = fragment_pool[i]->GetZ();
245 G4int A1 = fragment_pool[i]->GetA();
246 for(G4int j=0; j<nfrag; ++j) {
247 G4int Z2 = fragment_pool[j]->GetZ();
248 G4int A2 = fragment_pool[j]->GetA();
249 G4int Z = Z1 + Z2;
250 G4int A = A1 + A2;
251 if(Z < maxZ && A < maxA) {
252 if(IsAvailable(Z, A)){
253 std::vector<const G4VFermiFragment*> newvec;
254 newvec.push_back(fragment_pool[i]);
255 newvec.push_back(fragment_pool[j]);
256 if(!IsExist(Z, A, newvec)) {
257 G4FermiConfiguration* conf = new G4FermiConfiguration(newvec);
258 list2[A].push_back(conf);
259 ++counter;
260 }
261 }
262 }
263 }
264 }
265 if(verbose > 0) {
266 G4cout << G4endl;
267 G4cout << "### Pairs of fragments: " << counter << G4endl;
268 for(G4int A=2; A<maxA; ++A) {
269 G4cout << " A= " << A<<G4endl;
270 for(size_t j=0; j<list2[A].size(); ++j) {
271 std::vector<const G4VFermiFragment*> vector = (list2[A])[j]->GetFragmentList();
272 G4int a1=vector[0]->GetA();
273 G4int z1=vector[0]->GetZ();
274 G4int a2=vector[1]->GetA();
275 G4int z2=vector[1]->GetZ();
276 G4cout << "("<<a1<<","<<z1<<")("<<a2<<","<<z2<<") % ";
277 }
278 G4cout<<G4endl;
279 G4cout<<"---------------------------------------------------------------------------------"
280 << G4endl;
281 }
282 }
283
284 // list of fragment triples ordered by A
285 tot += counter;
286 counter = 0;
287 for(G4int A1=2; A1<maxA; ++A1) {
288 size_t nz = list2[A1].size();
289 for(size_t idx=0; idx<nz; ++idx) {
290 G4FermiConfiguration* conf2 = (list2[A1])[idx];
291 G4int Z1 = conf2->GetZ();
292 std::vector<const G4VFermiFragment*> vec2 = conf2->GetFragmentList();
293 //G4int a1 = vec2[0]->GetA();
294 // G4int z1 = vec2[0]->GetZ();
295 //G4int a2 = vec2[1]->GetA();
296 //G4int z2 = vec2[1]->GetZ();
297 for(G4int j=0; j<nfrag; ++j) {
298 G4int Z2 = fragment_pool[j]->GetZ();
299 G4int A2 = fragment_pool[j]->GetA();
300 G4int Z = Z1 + Z2;
301 G4int A = A1 + A2;
302 if(Z < maxZ && A < maxA) {
303 //if(IsAvailable(Z, A) && IsAvailable(z1+Z2, a1+A2)
304 // && IsAvailable(z2+Z2, a2+A2)) {
305 std::vector<const G4VFermiFragment*> newvec;
306 newvec.push_back(vec2[0]);
307 newvec.push_back(vec2[1]);
308 newvec.push_back(fragment_pool[j]);
309 if(!IsExist(Z, A, newvec)) {
310 G4FermiConfiguration* conf3 = new G4FermiConfiguration(newvec);
311 list3[A].push_back(conf3);
312 ++counter;
313 //}
314 }
315 }
316 }
317 }
318 }
319 if(verbose > 0) {
320 G4cout << G4endl;
321 G4cout << "### Triples of fragments: " << counter << G4endl;
322 for(G4int A=3; A<maxA; ++A) {
323 G4cout << " A= " << A<<G4endl;
324 for(size_t j=0; j<list3[A].size(); ++j) {
325 std::vector<const G4VFermiFragment*> vector = (list3[A])[j]->GetFragmentList();
326 G4int a1=vector[0]->GetA();
327 G4int z1=vector[0]->GetZ();
328 G4int a2=vector[1]->GetA();
329 G4int z2=vector[1]->GetZ();
330 G4int a3=vector[2]->GetA();
331 G4int z3=vector[2]->GetZ();
332 G4cout << "("<<a1<<","<<z1<<")("<<a2<<","<<z2<<")("<<a3<<","<<z3<<") % ";
333 }
334 G4cout<<G4endl;
335 G4cout<<"---------------------------------------------------------------------------------"
336 << G4endl;
337 }
338 }
339
340 // list of fragment quartets (3 + 1) ordered by A
341 tot += counter;
342 counter = 0;
343 for(G4int A1=3; A1<maxA; ++A1) {
344 size_t nz = list3[A1].size();
345 for(size_t idx=0; idx<nz; ++idx) {
346 G4FermiConfiguration* conf3 = (list3[A1])[idx];
347 G4int Z1 = conf3->GetZ();
348 std::vector<const G4VFermiFragment*> vec3 = conf3->GetFragmentList();
349 //G4int a1 = vec3[0]->GetA();
350 //G4int z1 = vec3[0]->GetZ();
351 //G4int a2 = vec3[1]->GetA();
352 //G4int z2 = vec3[1]->GetZ();
353 //G4int a3 = vec3[2]->GetA();
354 //G4int z3 = vec3[2]->GetZ();
355 for(G4int j=0; j<nfrag; ++j) {
356 G4int Z2 = fragment_pool[j]->GetZ();
357 G4int A2 = fragment_pool[j]->GetA();
358 G4int Z = Z1 + Z2;
359 G4int A = A1 + A2;
360 if(Z < maxZ && A < maxA) {
361 //if(IsAvailable(Z, A) && IsAvailable(z1+Z2, a1+A2)
362 // && IsAvailable(z2+Z2, a2+A2) && IsAvailable(z3+Z2, a3+A2)) {
363 std::vector<const G4VFermiFragment*> newvec;
364 newvec.push_back(vec3[0]);
365 newvec.push_back(vec3[1]);
366 newvec.push_back(vec3[2]);
367 newvec.push_back(fragment_pool[j]);
368 if(!IsExist(Z, A, newvec)) {
369 G4FermiConfiguration* conf4 = new G4FermiConfiguration(newvec);
370 list4[A].push_back(conf4);
371 ++counter;
372 }
373 //}
374 }
375 }
376 }
377 }
378 // list of fragment quartets (2 + 2) ordered by A
379 for(G4int A1=2; A1<maxA; ++A1) {
380 size_t nz1 = list2[A1].size();
381 for(size_t id1=0; id1<nz1; ++id1) {
382 G4FermiConfiguration* conf1 = (list2[A1])[id1];
383 G4int Z1 = conf1->GetZ();
384 std::vector<const G4VFermiFragment*> vec1 = conf1->GetFragmentList();
385 //G4int a1 = vec1[0]->GetA();
386 //G4int z1 = vec1[0]->GetZ();
387 //G4int a2 = vec1[1]->GetA();
388 //G4int z2 = vec1[1]->GetZ();
389 for(G4int A2=2; A2<maxA; ++A2) {
390 size_t nz2 = list2[A2].size();
391 for(size_t id2=0; id2<nz2; ++id2) {
392 G4FermiConfiguration* conf2 = (list2[A2])[id2];
393 G4int Z2 = conf2->GetZ();
394 std::vector<const G4VFermiFragment*> vec2 = conf2->GetFragmentList();
395 //G4int a3 = vec2[0]->GetA();
396 //G4int z3 = vec2[0]->GetZ();
397 //G4int a4 = vec2[1]->GetA();
398 //G4int z4 = vec2[1]->GetZ();
399 G4int Z = Z1 + Z2;
400 G4int A = A1 + A2;
401 if(Z < maxZ && A < maxA) {
402 //if(IsAvailable(Z, A) && IsAvailable(z1+z3, a1+a3)
403 // && IsAvailable(z1+z4, a1+a4) && IsAvailable(z2+z3, a2+a3)
404 // && IsAvailable(z2+z4, a2+a4) && IsAvailable(Z-z1, A-a1)
405 // && IsAvailable(Z-z2, A-a2) && IsAvailable(Z-z3, A-a3)) {
406 std::vector<const G4VFermiFragment*> newvec;
407 newvec.push_back(vec1[0]);
408 newvec.push_back(vec1[1]);
409 newvec.push_back(vec2[0]);
410 newvec.push_back(vec2[1]);
411 if(!IsExist(Z, A, newvec)) {
412 G4FermiConfiguration* conf4 = new G4FermiConfiguration(newvec);
413 list4[A].push_back(conf4);
414 ++counter;
415 }
416 //}
417 }
418 }
419 }
420 }
421 }
422 if(verbose > 0) {
423 tot += counter;
424 G4cout << G4endl;
425 G4cout << "### Quartets of fragments: " << counter << G4endl;
426 for(G4int A=4; A<maxA; ++A) {
427 G4cout << " A= " << A<<G4endl;
428 for(size_t j=0; j<list4[A].size(); ++j) {
429 std::vector<const G4VFermiFragment*> vector = (list4[A])[j]->GetFragmentList();
430 G4int a1=vector[0]->GetA();
431 G4int z1=vector[0]->GetZ();
432 G4int a2=vector[1]->GetA();
433 G4int z2=vector[1]->GetZ();
434 G4int a3=vector[2]->GetA();
435 G4int z3=vector[2]->GetZ();
436 G4int a4=vector[3]->GetA();
437 G4int z4=vector[3]->GetZ();
438
439 G4cout << "("<<a1<<","<<z1<<")("<<a2<<","<<z2<<")("<<a3<<","<<z3<<")("<<a4<<","<<z4<<") % ";
440 }
441 G4cout<<G4endl;
442 G4cout<<"---------------------------------------------------------------------------------"
443 << G4endl;
444 }
445 G4cout << "Total number: " << tot << G4endl;
446 }
447}
448
449const std::vector<G4FermiConfiguration*>*
450G4FermiFragmentsPool::GetConfigurationList(G4int Z, G4int A, G4double mass)
451{
452 //JMQ 040511 for printing the total number of configurations for a given A
453 G4int nconf=0;
454
455 std::vector<G4FermiConfiguration*>* v = new std::vector<G4FermiConfiguration*>;
456 if(Z >= maxZ || A >= maxA) { return v; }
457
458 //G4cout << "G4FermiFragmentsPool::GetConfigurationList:"
459 // << " Z= " << Z << " A= " << A << " Mass(GeV)= " << mass/GeV<< G4endl;
460
461 // look into pair list
462 size_t nz = list2[A].size();
463 if(0 < nz) {
464 for(size_t j=0; j<nz; ++j) {
465 G4FermiConfiguration* conf = (list2[A])[j];
466 if(Z == conf->GetZ() && mass >= conf->GetMass()) {
467 v->push_back(conf);
468 ++nconf;
469 }
470 //if(Z == conf->GetZ()) {
471 //G4cout << "Pair dM(MeV)= " << mass - conf->GetMass() << G4endl; }
472 }
473 }
474 // look into triple list
475 nz = list3[A].size();
476 if(0 < nz) {
477 for(size_t j=0; j<nz; ++j) {
478 G4FermiConfiguration* conf = (list3[A])[j];
479 if(Z == conf->GetZ() && mass >= conf->GetMass()) {
480 v->push_back(conf);
481 ++nconf;
482 }
483 //if(Z == conf->GetZ()) {
484 //G4cout << "Triple dM(MeV)= " << mass - conf->GetMass() << G4endl; }
485 }
486 }
487 // look into quartet list
488 nz = list4[A].size();
489 if(0 < nz) {
490 for(size_t j=0; j<nz; ++j) {
491 G4FermiConfiguration* conf = (list4[A])[j];
492 if(Z == conf->GetZ() && mass >= conf->GetMass()) {
493 v->push_back(conf);
494 ++nconf;
495 }
496 //if(Z == conf->GetZ()) {
497 // G4cout << "Quartet dM(MeV)= " << mass - conf->GetMass() << G4endl; }
498 }
499 }
500 // return if vector not empty
501 if(0 < v->size()) {
502 if(verbose > 0) {
503 G4double ExEn= mass - G4NucleiProperties::GetNuclearMass(A,Z);
504 G4cout<<"Total number of configurations = "<<nconf<<" for A= "
505 <<A<<" Z= "<<Z<<" E*= "<< ExEn<<" MeV"<<G4endl;
506 size_t size_vector_conf = v->size();
507 for(size_t jc=0; jc<size_vector_conf; ++jc) {
508 std::vector<const G4VFermiFragment*> v_frag = (*v)[jc]->GetFragmentList();
509 size_t size_vector_fragments = v_frag.size();
510 G4cout<<size_vector_fragments<<"-body configuration "<<jc+1<<": ";
511 for(size_t jf=0;jf<size_vector_fragments;++jf){
512 G4int af= v_frag[jf]->GetA();
513 G4int zf= v_frag[jf]->GetZ();
514 G4double ex=v_frag[jf]->GetExcitationEnergy();
515 G4cout<<"(a="<<af<<", z="<<zf<<", ex="<<ex<<") ";
516 }
517 G4cout<<G4endl;
518 G4cout<<"-----------------------------------------------------"<<G4endl;
519 }
520 }
521 return v;
522 }
523
524 // search in the pool and if found then return vector with one element
525 nz = list1[A].size();
526 G4FermiConfiguration* conf1 = 0;
527 if(0 < nz) {
528 for(size_t j=0; j<nz; ++j) {
529 G4FermiConfiguration* conf = (list1[A])[j];
530 //if(Z == conf->GetZ()) {
531 // G4cout << "Single dM(MeV)= " << mass - conf->GetMass() << G4endl; }
532
533 if(Z == conf->GetZ() && mass >= conf->GetMass()) {
534 if(!(conf->GetFragmentList())[0]->IsStable()) {
535 ++nconf;
536 v->push_back(conf);
537 if(verbose > 0) {
538 G4double ExEn= mass -G4NucleiProperties::GetNuclearMass(A,Z);
539 G4cout<<"Total number of configurations = "<<nconf<<" for A= "
540 <<A<<" Z= "<<Z<<" E*= "<< ExEn<<" MeV"<<G4endl;
541 size_t size_vector_conf=v->size();
542 for(size_t jc=0; jc<size_vector_conf; ++jc) {
543 std::vector<const G4VFermiFragment*> v_frag = (*v)[jc]->GetFragmentList();
544 size_t size_vector_fragments=v_frag.size();
545 G4cout<<"1 Fragment configuration "<<jc+1<<": ";
546 for(size_t jf=0;jf<size_vector_fragments;++jf){
547 G4int af= v_frag[jf]->GetA();
548 G4int zf= v_frag[jf]->GetZ();
549 G4double ex=v_frag[jf]->GetExcitationEnergy();
550 G4cout<<"(a="<<af<<", z="<<zf<<", ex="<<ex<<") ";
551 }
552 G4cout<<G4endl;
553 G4cout<<"-----------------------------------------------------"<<G4endl;
554 }
555 }
556 return v;
557 } else {
558 conf1 = conf;
559 break;
560 }
561 }
562 }
563 }
564
565 // search in the list of exotic configurations
566 nz = listextra.size();
567 if(0 < nz) {
568 for(size_t j=0; j<nz; ++j) {
569 G4FermiConfiguration* conf = listextra[j];
570 if(Z == conf->GetZ() && A == conf->GetA() &&
571 mass >= conf->GetMass()) {
572 ++nconf;
573 v->push_back(conf);
574 if(verbose > 0) {
575 G4double ExEn= mass -G4NucleiProperties::GetNuclearMass(A,Z);
576 G4cout<<"Total number of configurations = "<<nconf<<" for A= "
577 <<A<<" Z= "<<Z<<" E*= "<< ExEn<<" MeV"<<G4endl;
578 size_t size_vector_conf=v->size();
579 for(size_t jc=0; jc<size_vector_conf; ++jc) {
580 std::vector<const G4VFermiFragment*> v_frag = (*v)[jc]->GetFragmentList();
581 size_t size_vector_fragments=v_frag.size();
582 G4cout<<"Found exotic configuration -> configuration "<<jc+1<<": ";
583 for(size_t jf=0;jf<size_vector_fragments;++jf){
584 G4int af= v_frag[jf]->GetA();
585 G4int zf= v_frag[jf]->GetZ();
586 G4double ex=v_frag[jf]->GetExcitationEnergy();
587 G4cout<<"(a="<<af<<", z="<<zf<<", ex="<<ex<<") ";
588 }
589 G4cout<<G4endl;
590 G4cout<<"-----------------------------------------------------"<<G4endl;
591 }
592 }
593 return v;
594 }
595 }
596 }
597 //G4cout << "Explore dM(MeV)= "
598 // << mass - Z*proton_mass_c2 - (A-Z)*neutron_mass_c2 << G4endl;
599
600 // add new exotic configuration
601 if(mass > Z*proton_mass_c2 + (A-Z)*neutron_mass_c2) {
602 std::vector<const G4VFermiFragment*> newvec;
603 G4int idx = 1;
604 for(G4int i=0; i<A; ++i) {
605 if(i == Z) { idx = 0; }
606 newvec.push_back(fragment_pool[idx]);
607 }
608 G4FermiConfiguration* conf = new G4FermiConfiguration(newvec);
609 listextra.push_back(conf);
610 v->push_back(conf);
611 ++nconf;
612 if(verbose > 0) {
613 G4cout<<"Total number of configurations = "<<nconf<<G4cout;
614 G4double ExEn= mass -G4NucleiProperties::GetNuclearMass(A,Z);
615 G4cout<<"Total number of configurations = "<<nconf<<" for A= "
616 <<A<<" Z= "<<Z<<" E*= "<< ExEn<<" MeV"<<G4endl;
617 size_t size_vector_conf=v->size();
618 for(size_t jc=0; jc<size_vector_conf; ++jc) {
619 std::vector<const G4VFermiFragment*> v_frag = (*v)[jc]->GetFragmentList();
620 size_t size_vector_fragments=v_frag.size();
621 G4cout<<"New exotic configuration -> configuration "<<jc+1<<": ";
622 for(size_t jf=0;jf<size_vector_fragments;++jf){
623 G4int af= v_frag[jf]->GetA();
624 G4int zf= v_frag[jf]->GetZ();
625 G4double ex=v_frag[jf]->GetExcitationEnergy();
626 G4cout<<"(a="<<af<<", z="<<zf<<", ex="<<ex<<") ";
627 }
628 G4cout<<G4endl;
629 G4cout<<"-----------------------------------------------------"<<G4endl;
630 }
631 }
632 return v;
633 }
634
635 // only photon evaporation is possible
636 if(conf1) {
637 v->push_back(conf1);
638 ++nconf;
639 if(verbose > 0) {
640 G4cout<<"Total number of configurations = "<<nconf<<G4endl;
641 G4double ExEn= mass -G4NucleiProperties::GetNuclearMass(A,Z);
642 G4cout<<"Total number of configurations = "<<nconf<<" for A= "
643 <<A<<" Z= "<<Z<<" E*= "<< ExEn<<" MeV"<<G4endl;
644 size_t size_vector_conf=v->size();
645 for(size_t jc=0; jc<size_vector_conf; ++jc) {
646 std::vector<const G4VFermiFragment*> v_frag = (*v)[jc]->GetFragmentList();
647 size_t size_vector_fragments=v_frag.size();
648 G4cout<<"Only evaporation is possible -> configuration "<<jc+1<<": ";
649 for(size_t jf=0;jf<size_vector_fragments;++jf){
650 G4int af= v_frag[jf]->GetA();
651 G4int zf= v_frag[jf]->GetZ();
652 G4double ex=v_frag[jf]->GetExcitationEnergy();
653 G4cout<<"(a="<<af<<", z="<<zf<<", ex="<<ex<<") ";
654 }
655 G4cout<<G4endl;
656 G4cout<<"-----------------------------------------------------"<<G4endl;
657 }
658 }
659 return v;
660 }
661
662 //failer
663 if(verbose > 0) {
664 G4cout << "G4FermiFragmentsPool::GetConfigurationList: WARNING: not "
665 << "able decay fragment Z= " << Z << " A= " << A
666 << " Mass(GeV)= " << mass/GeV<< G4endl;
667 }
668 return v;
669}
670
671G4bool
672G4FermiFragmentsPool::IsExist(G4int Z, G4int A,
673 std::vector<const G4VFermiFragment*>& newconf)
674{
675 size_t nn = newconf.size();
676 G4double mass = 0.0;
677 for(size_t i=0; i<nn; ++i) { mass += newconf[i]->GetTotalEnergy(); }
678 // look into pair list
679 if(2 == nn) {
680 size_t nz = list2[A].size();
681 if(0 < nz) {
682 for(size_t j=0; j<nz; ++j) {
683 G4FermiConfiguration* conf = (list2[A])[j];
684 if(Z == conf->GetZ() && A == conf->GetA() &&
685 std::fabs(mass - conf->GetMass()) < keV) {return true; }
686 }
687 }
688 return false;
689 }
690 // look into triple list
691 if(3 == nn) {
692 size_t nz = list3[A].size();
693 if(0 < nz) {
694 for(size_t j=0; j<nz; ++j) {
695 G4FermiConfiguration* conf = (list3[A])[j];
696 if(Z == conf->GetZ() && A == conf->GetA() &&
697 std::fabs(mass - conf->GetMass()) < keV) { return true; }
698 }
699 }
700 return false;
701 }
702 // look into quartet list
703 if(4 == nn) {
704 size_t nz = list4[A].size();
705 if(0 < nz) {
706 for(size_t j=0; j<nz; ++j) {
707 G4FermiConfiguration* conf = (list4[A])[j];
708 if(Z == conf->GetZ() && A == conf->GetA() &&
709 std::fabs(mass - conf->GetMass()) < keV) { return true; }
710 }
711 }
712 return false;
713 }
714 return false;
715}
716
717const G4VFermiFragment*
718G4FermiFragmentsPool::GetFragment(G4int Z, G4int A)
719{
720 const G4VFermiFragment* f = 0;
721 if(Z >= maxZ || A >= maxA) { return f; }
722 size_t nz = list1[A].size();
723 if(0 < nz) {
724 for(size_t j=0; j<nz; ++j) {
725 G4FermiConfiguration* conf = (list1[A])[j];
726 if(Z == conf->GetZ()) { return (conf->GetFragmentList())[0]; }
727 }
728 }
729 return f;
730}
G4double
double G4double
Definition: G4Types.hh:64
G4int
int G4int
Definition: G4Types.hh:66
G4bool
bool G4bool
Definition: G4Types.hh:67
G4endl
#define G4endl
Definition: G4ios.hh:52
G4cout
G4DLLIMPORT std::ostream G4cout
G4FermiConfiguration::GetFragmentList
const std::vector< const G4VFermiFragment * > & GetFragmentList()
Definition: G4FermiConfiguration.hh:87
G4FermiFragmentsPool::GetConfigurationList
const std::vector< G4FermiConfiguration * > * GetConfigurationList(G4int Z, G4int A, G4double mass)
Definition: G4FermiFragmentsPool.cc:450
G4FermiFragmentsPool::Instance
static G4FermiFragmentsPool * Instance()
Definition: G4FermiFragmentsPool.cc:57
G4FermiFragmentsPool::GetFragment
const G4VFermiFragment * GetFragment(G4int Z, G4int A)
Definition: G4FermiFragmentsPool.cc:718
G4NucleiProperties::GetNuclearMass
static G4double GetNuclearMass(const G4double A, const G4double Z)
Definition: G4NucleiProperties.cc:53