| Geant4 Cross Reference |
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 // G4VUserPhysicsList implementation
27 //
28 // Original author: H.Kurashige (Kobe University), 9 January 1998
29 // --------------------------------------------------------------------
30
31 #include "G4VUserPhysicsList.hh"
32
33 #include "G4Material.hh"
34 #include "G4MaterialCutsCouple.hh"
35 #include "G4ParticleDefinition.hh"
36 #include "G4ParticleTable.hh"
37 #include "G4PhysicsListHelper.hh"
38 #include "G4ProcessManager.hh"
39 #include "G4ProductionCuts.hh"
40 #include "G4ProductionCutsTable.hh"
41 #include "G4Region.hh"
42 #include "G4RegionStore.hh"
43 #include "G4SystemOfUnits.hh"
44 #include "G4UImanager.hh"
45 #include "G4UnitsTable.hh"
46 #include "G4UserPhysicsListMessenger.hh"
47 #include "G4VTrackingManager.hh"
48 #include "G4ios.hh"
49 #include "globals.hh"
50
51 #include <fstream>
52 #include <iomanip>
53 #include <unordered_set>
54
55 // This static member is thread local. For each thread, it holds the array
56 // size of G4VUPLData instances.
57 //
58 #define G4MT_theMessenger ((this->subInstanceManager.offset[this->g4vuplInstanceID])._theMessenger)
59 #define G4MT_thePLHelper ((this->subInstanceManager.offset[this->g4vuplInstanceID])._thePLHelper)
60 #define fIsPhysicsTableBuilt \
61 ((this->subInstanceManager.offset[this->g4vuplInstanceID])._fIsPhysicsTableBuilt)
62 #define fDisplayThreshold \
63 ((this->subInstanceManager.offset[this->g4vuplInstanceID])._fDisplayThreshold)
64 #define theParticleIterator \
65 ((this->subInstanceManager.offset[this->g4vuplInstanceID])._theParticleIterator)
66
67 // This field helps to use the class G4VUPLManager
68 //
69 G4VUPLManager G4VUserPhysicsList::subInstanceManager;
70
71 // --------------------------------------------------------------------
72 void G4VUPLData::initialize()
73 {
74 _theParticleIterator = G4ParticleTable::GetParticleTable()->GetIterator();
75 _theMessenger = nullptr;
76 _thePLHelper = G4PhysicsListHelper::GetPhysicsListHelper();
77 _fIsPhysicsTableBuilt = false;
78 _fDisplayThreshold = 0;
79 }
80
81 // --------------------------------------------------------------------
82 G4VUserPhysicsList::G4VUserPhysicsList()
83 {
84 g4vuplInstanceID = subInstanceManager.CreateSubInstance(); // AND
85 // default cut value (1.0mm)
86 defaultCutValue = 1.0 * mm;
87
88 // pointer to the particle table
89 theParticleTable = G4ParticleTable::GetParticleTable();
90
91 // pointer to the cuts table
92 fCutsTable = G4ProductionCutsTable::GetProductionCutsTable();
93
94 // set energy range for SetCut calcuration
95 fCutsTable->SetEnergyRange(0.99 * keV, 100 * TeV);
96
97 // UI Messenger
98 G4MT_theMessenger = new G4UserPhysicsListMessenger(this); // AND
99
100 // PhysicsListHelper
101 G4MT_thePLHelper->SetVerboseLevel(verboseLevel); // AND
102
103 fIsPhysicsTableBuilt = false;
104 fDisplayThreshold = 0;
105 }
106
107 // --------------------------------------------------------------------
108 void G4VUserPhysicsList::InitializeWorker()
109 {
110 // Remember messengers are per-thread, so this needs to be done by each
111 // worker and due to the presence of "this" cannot be done in
112 // G4VUPLData::initialize()
113 G4MT_theMessenger = new G4UserPhysicsListMessenger(this);
114 }
115
116 // --------------------------------------------------------------------
117 void G4VUserPhysicsList::TerminateWorker()
118 {
119 RemoveProcessManager();
120 RemoveTrackingManager();
121 delete G4MT_theMessenger;
122 G4MT_theMessenger = nullptr;
123 }
124
125 // --------------------------------------------------------------------
126 G4VUserPhysicsList::~G4VUserPhysicsList()
127 {
128 delete G4MT_theMessenger;
129 G4MT_theMessenger = nullptr;
130
131 RemoveProcessManager();
132 RemoveTrackingManager();
133
134 // invoke DeleteAllParticle
135 theParticleTable->DeleteAllParticles();
136 }
137
138 // --------------------------------------------------------------------
139 G4VUserPhysicsList::G4VUserPhysicsList(const G4VUserPhysicsList& right)
140 : verboseLevel(right.verboseLevel),
141 defaultCutValue(right.defaultCutValue),
142 isSetDefaultCutValue(right.isSetDefaultCutValue),
143 fRetrievePhysicsTable(right.fRetrievePhysicsTable),
144 fStoredInAscii(right.fStoredInAscii),
145 fIsCheckedForRetrievePhysicsTable(right.fIsCheckedForRetrievePhysicsTable),
146 fIsRestoredCutValues(right.fIsRestoredCutValues),
147 directoryPhysicsTable(right.directoryPhysicsTable),
148 fDisableCheckParticleList(right.fDisableCheckParticleList)
149 {
150 g4vuplInstanceID = subInstanceManager.CreateSubInstance();
151 // pointer to the particle table
152 theParticleTable = G4ParticleTable::GetParticleTable();
153 theParticleIterator = theParticleTable->GetIterator();
154 // pointer to the cuts table
155 fCutsTable = G4ProductionCutsTable::GetProductionCutsTable();
156
157 // UI Messenger
158 G4MT_theMessenger = new G4UserPhysicsListMessenger(this);
159
160 // PhysicsListHelper
161 G4MT_thePLHelper = G4PhysicsListHelper::GetPhysicsListHelper();
162 G4MT_thePLHelper->SetVerboseLevel(verboseLevel);
163
164 fIsPhysicsTableBuilt =
165 right.GetSubInstanceManager().offset[right.GetInstanceID()]._fIsPhysicsTableBuilt;
166 fDisplayThreshold =
167 right.GetSubInstanceManager().offset[right.GetInstanceID()]._fDisplayThreshold;
168 }
169
170 // --------------------------------------------------------------------
171 G4VUserPhysicsList& G4VUserPhysicsList::operator=(const G4VUserPhysicsList& right)
172 {
173 if (this != &right) {
174 verboseLevel = right.verboseLevel;
175 defaultCutValue = right.defaultCutValue;
176 isSetDefaultCutValue = right.isSetDefaultCutValue;
177 fRetrievePhysicsTable = right.fRetrievePhysicsTable;
178 fStoredInAscii = right.fStoredInAscii;
179 fIsCheckedForRetrievePhysicsTable = right.fIsCheckedForRetrievePhysicsTable;
180 fIsRestoredCutValues = right.fIsRestoredCutValues;
181 directoryPhysicsTable = right.directoryPhysicsTable;
182 fIsPhysicsTableBuilt =
183 right.GetSubInstanceManager().offset[right.GetInstanceID()]._fIsPhysicsTableBuilt;
184 fDisplayThreshold =
185 right.GetSubInstanceManager().offset[right.GetInstanceID()]._fDisplayThreshold;
186 fDisableCheckParticleList = right.fDisableCheckParticleList;
187 }
188 return *this;
189 }
190
191 // --------------------------------------------------------------------
192 void G4VUserPhysicsList::AddProcessManager(G4ParticleDefinition* newParticle, G4ProcessManager*)
193 {
194 if (newParticle == nullptr) return;
195 G4Exception("G4VUserPhysicsList::AddProcessManager", "Run0252", JustWarning,
196 "This method is obsolete");
197 }
198
199 // --------------------------------------------------------------------
200 void G4VUserPhysicsList::InitializeProcessManager()
201 {
202 // Request lock for particle table access. Some changes are inside
203 // this critical region.
204 #ifdef G4MULTITHREADED
205 G4MUTEXLOCK(&G4ParticleTable::particleTableMutex());
206 G4ParticleTable::lockCount()++;
207 #endif
208 G4ParticleDefinition* gion = G4ParticleTable::GetParticleTable()->GetGenericIon();
209
210 // loop over all particles in G4ParticleTable
211 theParticleIterator->reset();
212 while ((*theParticleIterator)()) {
213 G4ParticleDefinition* particle = theParticleIterator->value();
214 G4ProcessManager* pmanager = particle->GetProcessManager();
215
216 if (pmanager == nullptr) {
217 // create process manager if the particle does not have its own.
218 pmanager = new G4ProcessManager(particle);
219 particle->SetProcessManager(pmanager);
220 if (particle->GetMasterProcessManager() == nullptr)
221 particle->SetMasterProcessManager(pmanager);
222 #ifdef G4VERBOSE
223 if (verboseLevel > 2) {
224 G4cout << "G4VUserPhysicsList::InitializeProcessManager: creating "
225 "ProcessManager to "
226 << particle->GetParticleName() << G4endl;
227 }
228 #endif
229 }
230 }
231
232 if (gion != nullptr) {
233 G4ProcessManager* gionPM = gion->GetProcessManager();
234 // loop over all particles once again (this time, with all general ions)
235 theParticleIterator->reset(false);
236 while ((*theParticleIterator)()) {
237 G4ParticleDefinition* particle = theParticleIterator->value();
238 if (particle->IsGeneralIon()) {
239 particle->SetProcessManager(gionPM);
240 #ifdef G4VERBOSE
241 if (verboseLevel > 2) {
242 G4cout << "G4VUserPhysicsList::InitializeProcessManager: copying "
243 "ProcessManager to "
244 << particle->GetParticleName() << G4endl;
245 }
246 #endif
247 }
248 }
249 }
250
251 // release lock for particle table access
252 #ifdef G4MULTITHREADED
253 G4MUTEXUNLOCK(&G4ParticleTable::particleTableMutex());
254 #endif
255 }
256
257 // --------------------------------------------------------------------
258 void G4VUserPhysicsList::RemoveProcessManager()
259 {
260 // Request lock for particle table access. Some changes are inside
261 // this critical region.
262 #ifdef G4MULTITHREADED
263 G4MUTEXLOCK(&G4ParticleTable::particleTableMutex());
264 G4ParticleTable::lockCount()++;
265 #endif
266
267 // loop over all particles in G4ParticleTable
268 theParticleIterator->reset();
269 while ((*theParticleIterator)()) {
270 G4ParticleDefinition* particle = theParticleIterator->value();
271 if (particle->GetInstanceID() < G4ParticleDefinitionSubInstanceManager::slavetotalspace()) {
272 if (particle->GetParticleSubType() != "generic"
273 || particle->GetParticleName() == "GenericIon")
274 {
275 G4ProcessManager* pmanager = particle->GetProcessManager();
276
277 delete pmanager;
278 #ifdef G4VERBOSE
279 if (verboseLevel > 2) {
280 G4cout << "G4VUserPhysicsList::RemoveProcessManager: ";
281 G4cout << "remove ProcessManager from ";
282 G4cout << particle->GetParticleName() << G4endl;
283 }
284 #endif
285 }
286 particle->SetProcessManager(nullptr);
287 }
288 }
289
290 // release lock for particle table access
291 #ifdef G4MULTITHREADED
292 G4MUTEXUNLOCK(&G4ParticleTable::particleTableMutex());
293 #endif
294 }
295
296 // --------------------------------------------------------------------
297 void G4VUserPhysicsList::RemoveTrackingManager()
298 {
299 // One tracking manager may be registered for multiple particles, make sure
300 // to delete every object only once.
301 std::unordered_set<G4VTrackingManager*> trackingManagers;
302
303 // loop over all particles in G4ParticleTable
304 theParticleIterator->reset();
305 while ((*theParticleIterator)()) {
306 G4ParticleDefinition* particle = theParticleIterator->value();
307 if (auto* trackingManager = particle->GetTrackingManager()) {
308 #ifdef G4VERBOSE
309 if (verboseLevel > 2) {
310 G4cout << "G4VUserPhysicsList::RemoveTrackingManager: ";
311 G4cout << "remove TrackingManager from ";
312 G4cout << particle->GetParticleName() << G4endl;
313 }
314 #endif
315 trackingManagers.insert(trackingManager);
316 particle->SetTrackingManager(nullptr);
317 }
318 }
319
320 for (G4VTrackingManager* tm : trackingManagers) {
321 delete tm;
322 }
323 }
324
325 // --------------------------------------------------------------------
326 void G4VUserPhysicsList::SetCuts()
327 {
328 if (!isSetDefaultCutValue) {
329 SetDefaultCutValue(defaultCutValue);
330 }
331
332 #ifdef G4VERBOSE
333 if (verboseLevel > 1) {
334 G4cout << "G4VUserPhysicsList::SetCuts: " << G4endl;
335 G4cout << "Cut for gamma: " << GetCutValue("gamma") / mm << "[mm]" << G4endl;
336 G4cout << "Cut for e-: " << GetCutValue("e-") / mm << "[mm]" << G4endl;
337 G4cout << "Cut for e+: " << GetCutValue("e+") / mm << "[mm]" << G4endl;
338 G4cout << "Cut for proton: " << GetCutValue("proton") / mm << "[mm]" << G4endl;
339 }
340 #endif
341
342 // dump Cut values if verboseLevel==3
343 if (verboseLevel > 2) {
344 DumpCutValuesTable();
345 }
346 }
347
348 // --------------------------------------------------------------------
349 void G4VUserPhysicsList::SetDefaultCutValue(G4double value)
350 {
351 if (value < 0.0) {
352 #ifdef G4VERBOSE
353 if (verboseLevel > 0) {
354 G4cout << "G4VUserPhysicsList::SetDefaultCutValue: negative cut values"
355 << " :" << value / mm << "[mm]" << G4endl;
356 }
357 #endif
358 return;
359 }
360
361 defaultCutValue = value;
362 isSetDefaultCutValue = true;
363
364 // set cut values for gamma at first and for e- and e+
365 SetCutValue(defaultCutValue, "gamma");
366 SetCutValue(defaultCutValue, "e-");
367 SetCutValue(defaultCutValue, "e+");
368 SetCutValue(defaultCutValue, "proton");
369
370 #ifdef G4VERBOSE
371 if (verboseLevel > 1) {
372 G4cout << "G4VUserPhysicsList::SetDefaultCutValue:"
373 << "default cut value is changed to :" << defaultCutValue / mm << "[mm]" << G4endl;
374 }
375 #endif
376 }
377
378 // --------------------------------------------------------------------
379 G4double G4VUserPhysicsList::GetCutValue(const G4String& name) const
380 {
381 std::size_t nReg = (G4RegionStore::GetInstance())->size();
382 if (nReg == 0) {
383 #ifdef G4VERBOSE
384 if (verboseLevel > 0) {
385 G4cout << "G4VUserPhysicsList::GetCutValue "
386 << " : No Default Region " << G4endl;
387 }
388 #endif
389 G4Exception("G4VUserPhysicsList::GetCutValue", "Run0253", FatalException, "No Default Region");
390 return -1. * mm;
391 }
392 G4Region* region = G4RegionStore::GetInstance()->GetRegion("DefaultRegionForTheWorld", false);
393 return region->GetProductionCuts()->GetProductionCut(name);
394 }
395
396 // --------------------------------------------------------------------
397 void G4VUserPhysicsList::SetCutValue(G4double aCut, const G4String& name)
398 {
399 SetParticleCuts(aCut, name);
400 }
401
402 // --------------------------------------------------------------------
403 void G4VUserPhysicsList::SetCutValue(G4double aCut, const G4String& pname, const G4String& rname)
404 {
405 G4Region* region = G4RegionStore::GetInstance()->GetRegion(rname);
406 if (region != nullptr) {
407 // set cut value
408 SetParticleCuts(aCut, pname, region);
409 }
410 else {
411 #ifdef G4VERBOSE
412 if (verboseLevel > 0) {
413 G4cout << "G4VUserPhysicsList::SetCutValue "
414 << " : No Region of " << rname << G4endl;
415 }
416 #endif
417 }
418 }
419
420 // --------------------------------------------------------------------
421 void G4VUserPhysicsList::SetCutsWithDefault()
422 {
423 SetDefaultCutValue(defaultCutValue);
424 G4VUserPhysicsList::SetCuts();
425 }
426
427 // --------------------------------------------------------------------
428 void G4VUserPhysicsList::SetCutsForRegion(G4double aCut, const G4String& rname)
429 {
430 // set cut values for gamma at first and for e- and e+
431 SetCutValue(aCut, "gamma", rname);
432 SetCutValue(aCut, "e-", rname);
433 SetCutValue(aCut, "e+", rname);
434 SetCutValue(aCut, "proton", rname);
435 }
436
437 // --------------------------------------------------------------------
438 void G4VUserPhysicsList::SetParticleCuts(G4double cut, G4ParticleDefinition* particle,
439 G4Region* region)
440 {
441 SetParticleCuts(cut, particle->GetParticleName(), region);
442 }
443
444 // --------------------------------------------------------------------
445 void G4VUserPhysicsList::SetParticleCuts(G4double cut, const G4String& particleName,
446 G4Region* region)
447 {
448 if (cut < 0.0) {
449 #ifdef G4VERBOSE
450 if (verboseLevel > 0) {
451 G4cout << "G4VUserPhysicsList::SetParticleCuts: negative cut values"
452 << " :" << cut / mm << "[mm]"
453 << " for " << particleName << G4endl;
454 }
455 #endif
456 return;
457 }
458
459 G4Region* world_region =
460 G4RegionStore::GetInstance()->GetRegion("DefaultRegionForTheWorld", false);
461 if (region == nullptr) {
462 std::size_t nReg = G4RegionStore::GetInstance()->size();
463 if (nReg == 0) {
464 #ifdef G4VERBOSE
465 if (verboseLevel > 0) {
466 G4cout << "G4VUserPhysicsList::SetParticleCuts "
467 << " : No Default Region " << G4endl;
468 }
469 #endif
470 G4Exception("G4VUserPhysicsList::SetParticleCuts ", "Run0254", FatalException,
471 "No Default Region");
472 return;
473 }
474 region = world_region;
475 }
476
477 if (!isSetDefaultCutValue) {
478 SetDefaultCutValue(defaultCutValue);
479 }
480
481 G4ProductionCuts* pcuts = region->GetProductionCuts();
482 if (region != world_region
483 && pcuts == G4ProductionCutsTable::GetProductionCutsTable()->GetDefaultProductionCuts())
484 { // This region had no unique cuts yet but shares the default cuts.
485 // Need to create a new object before setting the value.
486 pcuts = new G4ProductionCuts(
487 *(G4ProductionCutsTable::GetProductionCutsTable()->GetDefaultProductionCuts()));
488 region->SetProductionCuts(pcuts);
489 }
490 pcuts->SetProductionCut(cut, particleName);
491 #ifdef G4VERBOSE
492 if (verboseLevel > 2) {
493 G4cout << "G4VUserPhysicsList::SetParticleCuts: "
494 << " :" << cut / mm << "[mm]"
495 << " for " << particleName << G4endl;
496 }
497 #endif
498 }
499
500 // --------------------------------------------------------------------
501 void G4VUserPhysicsList::BuildPhysicsTable()
502 {
503 // Prepare Physics table for all particles
504 theParticleIterator->reset();
505 while ((*theParticleIterator)()) {
506 G4ParticleDefinition* particle = theParticleIterator->value();
507 PreparePhysicsTable(particle);
508 }
509
510 // ask processes to prepare physics table
511 if (fRetrievePhysicsTable) {
512 fIsRestoredCutValues = fCutsTable->RetrieveCutsTable(directoryPhysicsTable, fStoredInAscii);
513 // check if retrieve Cut Table successfully
514 if (!fIsRestoredCutValues) {
515 #ifdef G4VERBOSE
516 if (verboseLevel > 0) {
517 G4cout << "G4VUserPhysicsList::BuildPhysicsTable"
518 << " Retrieve Cut Table failed !!" << G4endl;
519 }
520 #endif
521 G4Exception("G4VUserPhysicsList::BuildPhysicsTable", "Run0255", RunMustBeAborted,
522 "Fail to retrieve Production Cut Table");
523 }
524 else {
525 #ifdef G4VERBOSE
526 if (verboseLevel > 2) {
527 G4cout << "G4VUserPhysicsList::BuildPhysicsTable"
528 << " Retrieve Cut Table successfully " << G4endl;
529 }
530 #endif
531 }
532 }
533 else {
534 #ifdef G4VERBOSE
535 if (verboseLevel > 2) {
536 G4cout << "G4VUserPhysicsList::BuildPhysicsTable"
537 << " does not retrieve Cut Table but calculate " << G4endl;
538 }
539 #endif
540 }
541
542 // Sets a value to particle
543 // set cut values for gamma at first and for e- and e+
544 G4String particleName;
545 G4ParticleDefinition* GammaP = theParticleTable->FindParticle("gamma");
546 if (GammaP != nullptr) BuildPhysicsTable(GammaP);
547 G4ParticleDefinition* EMinusP = theParticleTable->FindParticle("e-");
548 if (EMinusP != nullptr) BuildPhysicsTable(EMinusP);
549 G4ParticleDefinition* EPlusP = theParticleTable->FindParticle("e+");
550 if (EPlusP != nullptr) BuildPhysicsTable(EPlusP);
551 G4ParticleDefinition* ProtonP = theParticleTable->FindParticle("proton");
552 if (ProtonP != nullptr) BuildPhysicsTable(ProtonP);
553
554 theParticleIterator->reset();
555 while ((*theParticleIterator)()) {
556 G4ParticleDefinition* particle = theParticleIterator->value();
557 if (particle != GammaP && particle != EMinusP && particle != EPlusP && particle != ProtonP) {
558 BuildPhysicsTable(particle);
559 }
560 }
561
562 // Set flag
563 fIsPhysicsTableBuilt = true;
564 }
565
566 // --------------------------------------------------------------------
567 void G4VUserPhysicsList::BuildPhysicsTable(G4ParticleDefinition* particle)
568 {
569 if (auto* trackingManager = particle->GetTrackingManager()) {
570 #ifdef G4VERBOSE
571 if (verboseLevel > 2) {
572 G4cout << "G4VUserPhysicsList::BuildPhysicsTable "
573 << "Calculate Physics Table for " << particle->GetParticleName()
574 << " via custom TrackingManager" << G4endl;
575 }
576 #endif
577 trackingManager->BuildPhysicsTable(*particle);
578 return;
579 }
580
581 // Change in order to share physics tables for two kind of process.
582
583 if (particle->GetMasterProcessManager() == nullptr) {
584 #ifdef G4VERBOSE
585 if (verboseLevel > 0) {
586 G4cout << "#### G4VUserPhysicsList::BuildPhysicsTable() - BuildPhysicsTable("
587 << particle->GetParticleName() << ") skipped..." << G4endl;
588 }
589 #endif
590 return;
591 }
592 if (fRetrievePhysicsTable) {
593 if (!fIsRestoredCutValues) {
594 // fail to retrieve cut tables
595 #ifdef G4VERBOSE
596 if (verboseLevel > 0) {
597 G4cout << "G4VUserPhysicsList::BuildPhysicsTable "
598 << "Physics table can not be retrieved and will be calculated " << G4endl;
599 }
600 #endif
601 fRetrievePhysicsTable = false;
602 }
603 else {
604 #ifdef G4VERBOSE
605 if (verboseLevel > 2) {
606 G4cout << "G4VUserPhysicsList::BuildPhysicsTable "
607 << " Retrieve Physics Table for " << particle->GetParticleName() << G4endl;
608 }
609 #endif
610 // Retrieve PhysicsTable from files for proccesses
611 RetrievePhysicsTable(particle, directoryPhysicsTable, fStoredInAscii);
612 }
613 }
614
615 #ifdef G4VERBOSE
616 if (verboseLevel > 2) {
617 G4cout << "G4VUserPhysicsList::BuildPhysicsTable "
618 << "Calculate Physics Table for " << particle->GetParticleName() << G4endl;
619 }
620 #endif
621 // Rebuild the physics tables for every process for this particle type
622 // if particle is not ShortLived
623 if (!particle->IsShortLived()) {
624 G4ProcessManager* pManager = particle->GetProcessManager();
625 if (pManager == nullptr) {
626 #ifdef G4VERBOSE
627 if (verboseLevel > 0) {
628 G4cout << "G4VUserPhysicsList::BuildPhysicsTable "
629 << " : No Process Manager for " << particle->GetParticleName() << G4endl;
630 G4cout << particle->GetParticleName() << " should be created in your PhysicsList" << G4endl;
631 }
632 #endif
633 G4Exception("G4VUserPhysicsList::BuildPhysicsTable", "Run0271", FatalException,
634 "No process manager");
635 return;
636 }
637
638 // Get processes from master thread;
639 G4ProcessManager* pManagerShadow = particle->GetMasterProcessManager();
640
641 G4ProcessVector* pVector = pManager->GetProcessList();
642 if (pVector == nullptr) {
643 #ifdef G4VERBOSE
644 if (verboseLevel > 0) {
645 G4cout << "G4VUserPhysicsList::BuildPhysicsTable "
646 << " : No Process Vector for " << particle->GetParticleName() << G4endl;
647 }
648 #endif
649 G4Exception("G4VUserPhysicsList::BuildPhysicsTable", "Run0272", FatalException,
650 "No process Vector");
651 return;
652 }
653 #ifdef G4VERBOSE
654 if (verboseLevel > 2) {
655 G4cout << "G4VUserPhysicsList::BuildPhysicsTable %%%%%% " << particle->GetParticleName()
656 << G4endl;
657 G4cout << " ProcessManager : " << pManager << " ProcessManagerShadow : " << pManagerShadow
658 << G4endl;
659 for (G4int iv1 = 0; iv1 < (G4int)pVector->size(); ++iv1) {
660 G4cout << " " << iv1 << " - " << (*pVector)[iv1]->GetProcessName() << G4endl;
661 }
662 G4cout << "--------------------------------------------------------------" << G4endl;
663 G4ProcessVector* pVectorShadow = pManagerShadow->GetProcessList();
664
665 for (G4int iv2 = 0; iv2 < (G4int)pVectorShadow->size(); ++iv2) {
666 G4cout << " " << iv2 << " - " << (*pVectorShadow)[iv2]->GetProcessName() << G4endl;
667 }
668 }
669 #endif
670 for (G4int j = 0; j < (G4int)pVector->size(); ++j) {
671 // Andrea July 16th 2013 : migration to new interface...
672 // Infer if we are in a worker thread or master thread
673 // Master thread is the one in which the process manager
674 // and process manager shadow pointers are the same
675 if (pManagerShadow == pManager) {
676 (*pVector)[j]->BuildPhysicsTable(*particle);
677 }
678 else {
679 (*pVector)[j]->BuildWorkerPhysicsTable(*particle);
680 }
681
682 } // End loop on processes vector
683 } // End if short-lived
684 }
685
686 // --------------------------------------------------------------------
687 void G4VUserPhysicsList::PreparePhysicsTable(G4ParticleDefinition* particle)
688 {
689 if (auto* trackingManager = particle->GetTrackingManager()) {
690 trackingManager->PreparePhysicsTable(*particle);
691 return;
692 }
693
694 if ((particle->GetMasterProcessManager()) == nullptr) {
695 return;
696 }
697 // Prepare the physics tables for every process for this particle type
698 // if particle is not ShortLived
699 if (!particle->IsShortLived()) {
700 G4ProcessManager* pManager = particle->GetProcessManager();
701 if (pManager == nullptr) {
702 #ifdef G4VERBOSE
703 if (verboseLevel > 0) {
704 G4cout << "G4VUserPhysicsList::PreparePhysicsTable "
705 << ": No Process Manager for " << particle->GetParticleName() << G4endl;
706 G4cout << particle->GetParticleName() << " should be created in your PhysicsList" << G4endl;
707 }
708 #endif
709 G4Exception("G4VUserPhysicsList::PreparePhysicsTable", "Run0273", FatalException,
710 "No process manager");
711 return;
712 }
713
714 // Get processes from master thread
715 G4ProcessManager* pManagerShadow = particle->GetMasterProcessManager();
716
717 G4ProcessVector* pVector = pManager->GetProcessList();
718 if (pVector == nullptr) {
719 #ifdef G4VERBOSE
720 if (verboseLevel > 0) {
721 G4cout << "G4VUserPhysicsList::PreparePhysicsTable "
722 << ": No Process Vector for " << particle->GetParticleName() << G4endl;
723 }
724 #endif
725 G4Exception("G4VUserPhysicsList::PreparePhysicsTable", "Run0274", FatalException,
726 "No process Vector");
727 return;
728 }
729 for (G4int j = 0; j < (G4int)pVector->size(); ++j) {
730 // Andrea July 16th 2013 : migration to new interface...
731 // Infer if we are in a worker thread or master thread
732 // Master thread is the one in which the process manager
733 // and process manager shadow pointers are the same
734 if (pManagerShadow == pManager) {
735 (*pVector)[j]->PreparePhysicsTable(*particle);
736 }
737 else {
738 (*pVector)[j]->PrepareWorkerPhysicsTable(*particle);
739 }
740 } // End loop on processes vector
741 } // End if pn ShortLived
742 }
743
744 // --------------------------------------------------------------------
745 void G4VUserPhysicsList::BuildIntegralPhysicsTable(G4VProcess* process,
746 G4ParticleDefinition* particle)
747 {
748 // TODO Should we change this function?
749 //*******************************************************************
750 // Temporary addition to make the integral schema of electromagnetic
751 // processes work.
752
753 if ((process->GetProcessName() == "Imsc") || (process->GetProcessName() == "IeIoni")
754 || (process->GetProcessName() == "IeBrems") || (process->GetProcessName() == "Iannihil")
755 || (process->GetProcessName() == "IhIoni") || (process->GetProcessName() == "IMuIoni")
756 || (process->GetProcessName() == "IMuBrems") || (process->GetProcessName() == "IMuPairProd"))
757 {
758 #ifdef G4VERBOSE
759 if (verboseLevel > 2) {
760 G4cout << "G4VUserPhysicsList::BuildIntegralPhysicsTable "
761 << " BuildPhysicsTable is invoked for " << process->GetProcessName() << "("
762 << particle->GetParticleName() << ")" << G4endl;
763 }
764 #endif
765 process->BuildPhysicsTable(*particle);
766 }
767 }
768
769 // --------------------------------------------------------------------
770 void G4VUserPhysicsList::DumpList() const
771 {
772 theParticleIterator->reset();
773 G4int idx = 0;
774 while ((*theParticleIterator)()) {
775 G4ParticleDefinition* particle = theParticleIterator->value();
776 G4cout << particle->GetParticleName();
777 if ((idx++ % 4) == 3) {
778 G4cout << G4endl;
779 }
780 else {
781 G4cout << ", ";
782 }
783 }
784 G4cout << G4endl;
785 }
786
787 // --------------------------------------------------------------------
788 void G4VUserPhysicsList::DumpCutValuesTable(G4int flag)
789 {
790 fDisplayThreshold = flag;
791 }
792
793 // --------------------------------------------------------------------
794 void G4VUserPhysicsList::DumpCutValuesTableIfRequested()
795 {
796 if (fDisplayThreshold == 0) return;
797 G4ProductionCutsTable::GetProductionCutsTable()->DumpCouples();
798 fDisplayThreshold = 0;
799 }
800
801 // --------------------------------------------------------------------
802 G4bool G4VUserPhysicsList::StorePhysicsTable(const G4String& directory)
803 {
804 G4bool ascii = fStoredInAscii;
805 G4String dir = directory;
806 if (dir.empty())
807 dir = directoryPhysicsTable;
808 else
809 directoryPhysicsTable = dir;
810
811 // store CutsTable info
812 if (!fCutsTable->StoreCutsTable(dir, ascii)) {
813 G4Exception("G4VUserPhysicsList::StorePhysicsTable", "Run0281", JustWarning,
814 "Fail to store Cut Table");
815 return false;
816 }
817 #ifdef G4VERBOSE
818 if (verboseLevel > 2) {
819 G4cout << "G4VUserPhysicsList::StorePhysicsTable "
820 << " Store material and cut values successfully" << G4endl;
821 }
822 #endif
823
824 G4bool success = true;
825
826 // loop over all particles in G4ParticleTable
827 theParticleIterator->reset();
828 while ((*theParticleIterator)()) {
829 G4ParticleDefinition* particle = theParticleIterator->value();
830 // Store physics tables for every process for this particle type
831 G4ProcessVector* pVector = (particle->GetProcessManager())->GetProcessList();
832 for (G4int j = 0; j < (G4int)pVector->size(); ++j) {
833 if (!(*pVector)[j]->StorePhysicsTable(particle, dir, ascii)) {
834 G4String comment = "Fail to store physics table for ";
835 comment += (*pVector)[j]->GetProcessName();
836 comment += "(" + particle->GetParticleName() + ")";
837 G4Exception("G4VUserPhysicsList::StorePhysicsTable", "Run0282", JustWarning, comment);
838 success = false;
839 }
840 }
841 // end loop over processes
842 }
843 // end loop over particles
844 return success;
845 }
846
847 // --------------------------------------------------------------------
848 void G4VUserPhysicsList::SetPhysicsTableRetrieved(const G4String& directory)
849 {
850 fRetrievePhysicsTable = true;
851 if (!directory.empty()) {
852 directoryPhysicsTable = directory;
853 }
854 fIsCheckedForRetrievePhysicsTable = false;
855 fIsRestoredCutValues = false;
856 }
857
858 // --------------------------------------------------------------------
859 void G4VUserPhysicsList::RetrievePhysicsTable(G4ParticleDefinition* particle,
860 const G4String& directory, G4bool ascii)
861 {
862 G4bool success[100];
863 // Retrieve physics tables for every process for this particle type
864 G4ProcessVector* pVector = (particle->GetProcessManager())->GetProcessList();
865 for (G4int j = 0; j < (G4int)pVector->size(); ++j) {
866 success[j] = (*pVector)[j]->RetrievePhysicsTable(particle, directory, ascii);
867
868 if (!success[j]) {
869 #ifdef G4VERBOSE
870 if (verboseLevel > 2) {
871 G4cout << "G4VUserPhysicsList::RetrievePhysicsTable "
872 << " Fail to retrieve Physics Table for " << (*pVector)[j]->GetProcessName()
873 << G4endl;
874 G4cout << "Calculate Physics Table for " << particle->GetParticleName() << G4endl;
875 }
876 #endif
877 (*pVector)[j]->BuildPhysicsTable(*particle);
878 }
879 }
880 for (G4int j = 0; j < (G4int)pVector->size(); ++j) {
881 // temporary addition to make the integral schema
882 if (!success[j]) BuildIntegralPhysicsTable((*pVector)[j], particle);
883 }
884 }
885
886 // --------------------------------------------------------------------
887 void G4VUserPhysicsList::SetApplyCuts(G4bool value, const G4String& name)
888 {
889 #ifdef G4VERBOSE
890 if (verboseLevel > 2) {
891 G4cout << "G4VUserPhysicsList::SetApplyCuts for " << name << G4endl;
892 }
893 #endif
894 if (name == "all") {
895 theParticleTable->FindParticle("gamma")->SetApplyCutsFlag(value);
896 theParticleTable->FindParticle("e-")->SetApplyCutsFlag(value);
897 theParticleTable->FindParticle("e+")->SetApplyCutsFlag(value);
898 theParticleTable->FindParticle("proton")->SetApplyCutsFlag(value);
899 }
900 else {
901 theParticleTable->FindParticle(name)->SetApplyCutsFlag(value);
902 }
903 }
904
905 // --------------------------------------------------------------------
906 G4bool G4VUserPhysicsList::GetApplyCuts(const G4String& name) const
907 {
908 return theParticleTable->FindParticle(name)->GetApplyCutsFlag();
909 }
910
911 // --------------------------------------------------------------------
912 void G4VUserPhysicsList::CheckParticleList()
913 {
914 if (!fDisableCheckParticleList) {
915 G4MT_thePLHelper->CheckParticleList();
916 }
917 }
918
919 // --------------------------------------------------------------------
920 void G4VUserPhysicsList::AddTransportation()
921 {
922 G4MT_thePLHelper->AddTransportation();
923 }
924
925 // --------------------------------------------------------------------
926 void G4VUserPhysicsList::UseCoupledTransportation(G4bool vl)
927 {
928 G4MT_thePLHelper->UseCoupledTransportation(vl);
929 }
930
931 // --------------------------------------------------------------------
932 G4bool G4VUserPhysicsList::RegisterProcess(G4VProcess* process, G4ParticleDefinition* particle)
933 {
934 return G4MT_thePLHelper->RegisterProcess(process, particle);
935 }
936
937 // --------------------------------------------------------------------
938 G4ParticleTable::G4PTblDicIterator* G4VUserPhysicsList::GetParticleIterator() const
939 {
940 return (subInstanceManager.offset[g4vuplInstanceID])._theParticleIterator;
941 }
942
943 // --------------------------------------------------------------------
944 void G4VUserPhysicsList::SetVerboseLevel(G4int value)
945 {
946 verboseLevel = value;
947 // set verboseLevel for G4ProductionCutsTable same as one for
948 // G4VUserPhysicsList:
949 fCutsTable->SetVerboseLevel(verboseLevel);
950
951 G4MT_thePLHelper->SetVerboseLevel(verboseLevel);
952
953 #ifdef G4VERBOSE
954 if (verboseLevel > 1) {
955 G4cout << "G4VUserPhysicsList::SetVerboseLevel :"
956 << " Verbose level is set to " << verboseLevel << G4endl;
957 }
958 #endif
959 }
960