Geant4 9.6.0
Toolkit for the simulation of the passage of particles through matter
Loading...
Searching...
No Matches
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