Geant4-11
Public Member Functions | Static Public Member Functions | Protected Member Functions | Protected Attributes | Static Protected Attributes | Private Types
G4VUserPhysicsList Class Referenceabstract

#include <G4VUserPhysicsList.hh>

Inheritance diagram for G4VUserPhysicsList:
ExN01PhysicsList G4ErrorPhysicsList G4VModularPhysicsList PhysicsListEMstd pyG4VUserPhysicsList::CB_G4VUserPhysicsList FTFP_BERT FTFP_BERT_ATL FTFP_BERT_HP FTFP_BERT_TRV FTFQGSP_BERT FTF_BIC LBE NuBeam PhysicsList QBBC QGSP_BERT QGSP_BERT_HP QGSP_BIC QGSP_BIC_AllHP QGSP_BIC_HP QGSP_FTFP_BERT QGS_BIC QPhysicsList Shielding pyG4VModularPhysicsList::CB_G4VModularPhysicsList

Public Member Functions

void AddProcessManager (G4ParticleDefinition *newParticle, G4ProcessManager *newManager=nullptr)
 
void BuildPhysicsTable ()
 
void BuildPhysicsTable (G4ParticleDefinition *)
 
void CheckParticleList ()
 
void Construct ()
 
virtual void ConstructParticle ()=0
 
virtual void ConstructProcess ()=0
 
void DisableCheckParticleList ()
 
void DumpCutValuesTable (G4int flag=1)
 
void DumpCutValuesTableIfRequested ()
 
void DumpList () const
 
 G4VUserPhysicsList ()
 
 G4VUserPhysicsList (const G4VUserPhysicsList &)
 
G4bool GetApplyCuts (const G4String &name) const
 
G4double GetCutValue (const G4String &pname) const
 
G4double GetDefaultCutValue () const
 
G4int GetInstanceID () const
 
const G4StringGetPhysicsTableDirectory () const
 
G4int GetVerboseLevel () const
 
virtual void InitializeWorker ()
 
G4bool IsPhysicsTableRetrieved () const
 
G4bool IsStoredInAscii () const
 
G4VUserPhysicsListoperator= (const G4VUserPhysicsList &)
 
void PreparePhysicsTable (G4ParticleDefinition *)
 
void RemoveProcessManager ()
 
void RemoveTrackingManager ()
 
void ResetPhysicsTableRetrieved ()
 
void ResetStoredInAscii ()
 
void SetApplyCuts (G4bool value, const G4String &name)
 
virtual void SetCuts ()
 
void SetCutsForRegion (G4double aCut, const G4String &rname)
 
void SetCutsWithDefault ()
 
void SetCutValue (G4double aCut, const G4String &pname)
 
void SetCutValue (G4double aCut, const G4String &pname, const G4String &rname)
 
void SetDefaultCutValue (G4double newCutValue)
 
void SetParticleCuts (G4double cut, const G4String &particleName, G4Region *region=nullptr)
 
void SetParticleCuts (G4double cut, G4ParticleDefinition *particle, G4Region *region=nullptr)
 
void SetPhysicsTableRetrieved (const G4String &directory="")
 
void SetStoredInAscii ()
 
void SetVerboseLevel (G4int value)
 
G4bool StorePhysicsTable (const G4String &directory=".")
 
virtual void TerminateWorker ()
 
void UseCoupledTransportation (G4bool vl=true)
 
virtual ~G4VUserPhysicsList ()
 

Static Public Member Functions

static const G4VUPLManagerGetSubInstanceManager ()
 

Protected Member Functions

void AddTransportation ()
 
void BuildIntegralPhysicsTable (G4VProcess *, G4ParticleDefinition *)
 
G4ParticleTable::G4PTblDicIteratorGetParticleIterator () const
 
void InitializeProcessManager ()
 
G4bool RegisterProcess (G4VProcess *process, G4ParticleDefinition *particle)
 
virtual void RetrievePhysicsTable (G4ParticleDefinition *, const G4String &directory, G4bool ascii=false)
 

Protected Attributes

G4double defaultCutValue = 1.0
 
G4String directoryPhysicsTable = "."
 
G4ProductionCutsTablefCutsTable = nullptr
 
G4bool fDisableCheckParticleList = false
 
G4bool fIsCheckedForRetrievePhysicsTable = false
 
G4bool fIsRestoredCutValues = false
 
G4bool fRetrievePhysicsTable = false
 
G4bool fStoredInAscii = true
 
G4int g4vuplInstanceID = 0
 
G4bool isSetDefaultCutValue = false
 
G4ParticleTabletheParticleTable = nullptr
 
G4int verboseLevel = 1
 

Static Protected Attributes

static G4RUN_DLL G4VUPLManager subInstanceManager
 

Private Types

enum  { FixedStringLengthForStore = 32 }
 

Detailed Description

Definition at line 106 of file G4VUserPhysicsList.hh.

Member Enumeration Documentation

◆ anonymous enum

anonymous enum
private
Enumerator
FixedStringLengthForStore 

Definition at line 311 of file G4VUserPhysicsList.hh.

312 {
314 };

Constructor & Destructor Documentation

◆ G4VUserPhysicsList() [1/2]

G4VUserPhysicsList::G4VUserPhysicsList ( )

Definition at line 86 of file G4VUserPhysicsList.cc.

87{
89 // default cut value (1.0mm)
90 defaultCutValue = 1.0 * mm;
91
92 // pointer to the particle table
94 // theParticleIterator = theParticleTable->GetIterator();
95
96 // pointer to the cuts table
98
99 // set energy range for SetCut calcuration
100 fCutsTable->SetEnergyRange(0.99 * keV, 100 * TeV);
101
102 // UI Messenger
103 // theMessenger = new G4UserPhysicsListMessenger(this);
105
106 // PhysicsListHelper
107 // thePLHelper = G4PhysicsListHelper::GetPhysicsListHelper();
108 // thePLHelper->SetVerboseLevel(verboseLevel);
109 // G4MT_thePLHelper = G4PhysicsListHelper::GetPhysicsListHelper(); //AND
110 G4MT_thePLHelper->SetVerboseLevel(verboseLevel); // AND
111
112 fIsPhysicsTableBuilt = false;
114}
static constexpr double mm
Definition: G4SIunits.hh:95
static constexpr double keV
Definition: G4SIunits.hh:202
static constexpr double TeV
Definition: G4SIunits.hh:204
#define G4MT_theMessenger
#define G4MT_thePLHelper
#define fIsPhysicsTableBuilt
#define fDisplayThreshold
static G4ParticleTable * GetParticleTable()
void SetEnergyRange(G4double lowedge, G4double highedge)
static G4ProductionCutsTable * GetProductionCutsTable()
G4int CreateSubInstance()
G4ProductionCutsTable * fCutsTable
G4ParticleTable * theParticleTable
static G4RUN_DLL G4VUPLManager subInstanceManager

References G4VUPLSplitter< T >::CreateSubInstance(), defaultCutValue, fCutsTable, fDisplayThreshold, fIsPhysicsTableBuilt, G4MT_theMessenger, G4MT_thePLHelper, g4vuplInstanceID, G4ParticleTable::GetParticleTable(), G4ProductionCutsTable::GetProductionCutsTable(), keV, mm, G4ProductionCutsTable::SetEnergyRange(), subInstanceManager, TeV, theParticleTable, and verboseLevel.

◆ ~G4VUserPhysicsList()

G4VUserPhysicsList::~G4VUserPhysicsList ( )
virtual

◆ G4VUserPhysicsList() [2/2]

G4VUserPhysicsList::G4VUserPhysicsList ( const G4VUserPhysicsList right)

Definition at line 148 of file G4VUserPhysicsList.cc.

158{
160 // pointer to the particle table
163 // pointer to the cuts table
165
166 // UI Messenger
168
169 // PhysicsListHelper
171 G4MT_thePLHelper->SetVerboseLevel(verboseLevel);
172
174 .offset[right.GetInstanceID()]
175 ._fIsPhysicsTableBuilt;
177 .offset[right.GetInstanceID()]
178 ._fDisplayThreshold;
179}
#define theParticleIterator
G4PTblDicIterator * GetIterator() const
static G4PhysicsListHelper * GetPhysicsListHelper()
G4RUN_DLL G4ThreadLocalStatic T * offset
G4bool fIsCheckedForRetrievePhysicsTable
G4int GetInstanceID() const
static const G4VUPLManager & GetSubInstanceManager()

References G4VUPLSplitter< T >::CreateSubInstance(), fCutsTable, fDisplayThreshold, fIsPhysicsTableBuilt, G4MT_theMessenger, G4MT_thePLHelper, g4vuplInstanceID, GetInstanceID(), G4ParticleTable::GetIterator(), G4ParticleTable::GetParticleTable(), G4PhysicsListHelper::GetPhysicsListHelper(), G4ProductionCutsTable::GetProductionCutsTable(), GetSubInstanceManager(), G4VUPLSplitter< T >::offset, subInstanceManager, theParticleIterator, theParticleTable, and verboseLevel.

Member Function Documentation

◆ AddProcessManager()

void G4VUserPhysicsList::AddProcessManager ( G4ParticleDefinition newParticle,
G4ProcessManager newManager = nullptr 
)

Definition at line 207 of file G4VUserPhysicsList.cc.

209{
210 if(newParticle == nullptr)
211 return;
212 G4Exception("G4VUserPhysicsList::AddProcessManager", "Run0252",
213 JustWarning, "This method is obsolete");
214}
@ JustWarning
void G4Exception(const char *originOfException, const char *exceptionCode, G4ExceptionSeverity severity, const char *description)
Definition: G4Exception.cc:35

References G4Exception(), and JustWarning.

Referenced by G4UserPhysicsListMessenger::SetNewValue().

◆ AddTransportation()

void G4VUserPhysicsList::AddTransportation ( )
protected

◆ BuildIntegralPhysicsTable()

void G4VUserPhysicsList::BuildIntegralPhysicsTable ( G4VProcess process,
G4ParticleDefinition particle 
)
protected

Definition at line 864 of file G4VUserPhysicsList.cc.

866{
867 // TODO Should we change this function?
868 //*******************************************************************
869 // Temporary addition to make the integral schema of electromagnetic
870 // processes work.
871
872 if((process->GetProcessName() == "Imsc") ||
873 (process->GetProcessName() == "IeIoni") ||
874 (process->GetProcessName() == "IeBrems") ||
875 (process->GetProcessName() == "Iannihil") ||
876 (process->GetProcessName() == "IhIoni") ||
877 (process->GetProcessName() == "IMuIoni") ||
878 (process->GetProcessName() == "IMuBrems") ||
879 (process->GetProcessName() == "IMuPairProd"))
880 {
881#ifdef G4VERBOSE
882 if(verboseLevel > 2)
883 {
884 G4cout << "G4VUserPhysicsList::BuildIntegralPhysicsTable "
885 << " BuildPhysicsTable is invoked for "
886 << process->GetProcessName() << "(" << particle->GetParticleName()
887 << ")" << G4endl;
888 }
889#endif
890 process->BuildPhysicsTable(*particle);
891 }
892}
#define G4endl
Definition: G4ios.hh:57
G4GLOB_DLL std::ostream G4cout
const G4String & GetParticleName() const
virtual void BuildPhysicsTable(const G4ParticleDefinition &)
Definition: G4VProcess.hh:187
const G4String & GetProcessName() const
Definition: G4VProcess.hh:382

References G4VProcess::BuildPhysicsTable(), G4cout, G4endl, G4ParticleDefinition::GetParticleName(), G4VProcess::GetProcessName(), and verboseLevel.

Referenced by RetrievePhysicsTable().

◆ BuildPhysicsTable() [1/2]

void G4VUserPhysicsList::BuildPhysicsTable ( )

Definition at line 562 of file G4VUserPhysicsList.cc.

563{
564 // Prepare Physics table for all particles
565 theParticleIterator->reset();
566 while((*theParticleIterator)())
567 {
568 G4ParticleDefinition* particle = theParticleIterator->value();
569 PreparePhysicsTable(particle);
570 }
571
572 // ask processes to prepare physics table
574 {
577 // check if retrieve Cut Table successfully
579 {
580#ifdef G4VERBOSE
581 if(verboseLevel > 0)
582 {
583 G4cout << "G4VUserPhysicsList::BuildPhysicsTable"
584 << " Retrieve Cut Table failed !!" << G4endl;
585 }
586#endif
587 G4Exception("G4VUserPhysicsList::BuildPhysicsTable", "Run0255",
588 RunMustBeAborted, "Fail to retrieve Production Cut Table");
589 }
590 else
591 {
592#ifdef G4VERBOSE
593 if(verboseLevel > 2)
594 {
595 G4cout << "G4VUserPhysicsList::BuildPhysicsTable"
596 << " Retrieve Cut Table successfully " << G4endl;
597 }
598#endif
599 }
600 }
601 else
602 {
603#ifdef G4VERBOSE
604 if(verboseLevel > 2)
605 {
606 G4cout << "G4VUserPhysicsList::BuildPhysicsTable"
607 << " does not retrieve Cut Table but calculate " << G4endl;
608 }
609#endif
610 }
611
612 // Sets a value to particle
613 // set cut values for gamma at first and for e- and e+
614 G4String particleName;
616 if(GammaP)
617 BuildPhysicsTable(GammaP);
619 if(EMinusP)
620 BuildPhysicsTable(EMinusP);
622 if(EPlusP)
623 BuildPhysicsTable(EPlusP);
625 if(ProtonP)
626 BuildPhysicsTable(ProtonP);
627
628 theParticleIterator->reset();
629 while((*theParticleIterator)())
630 {
631 G4ParticleDefinition* particle = theParticleIterator->value();
632 if(particle != GammaP && particle != EMinusP && particle != EPlusP &&
633 particle != ProtonP)
634 {
635 BuildPhysicsTable(particle);
636 }
637 }
638
639 // Set flag
641}
@ RunMustBeAborted
G4ParticleDefinition * FindParticle(G4int PDGEncoding)
G4bool RetrieveCutsTable(const G4String &directory, G4bool ascii=false)
void PreparePhysicsTable(G4ParticleDefinition *)

References BuildPhysicsTable(), directoryPhysicsTable, fCutsTable, G4ParticleTable::FindParticle(), fIsPhysicsTableBuilt, fIsRestoredCutValues, fRetrievePhysicsTable, fStoredInAscii, G4cout, G4endl, G4Exception(), PreparePhysicsTable(), G4ProductionCutsTable::RetrieveCutsTable(), RunMustBeAborted, theParticleIterator, theParticleTable, and verboseLevel.

Referenced by BuildPhysicsTable(), G4RunManagerKernel::BuildPhysicsTables(), and G4UserPhysicsListMessenger::SetNewValue().

◆ BuildPhysicsTable() [2/2]

void G4VUserPhysicsList::BuildPhysicsTable ( G4ParticleDefinition particle)

Definition at line 644 of file G4VUserPhysicsList.cc.

645{
646 if (auto *trackingManager = particle->GetTrackingManager())
647 {
648#ifdef G4VERBOSE
649 if(verboseLevel > 2)
650 {
651 G4cout << "G4VUserPhysicsList::BuildPhysicsTable "
652 << "Calculate Physics Table for " << particle->GetParticleName()
653 << " via custom TrackingManager" << G4endl;
654 }
655#endif
656 trackingManager->BuildPhysicsTable(*particle);
657 return;
658 }
659
660 // Change in order to share physics tables for two kind of process.
661
662 if(particle->GetMasterProcessManager() == nullptr)
663 {
664#ifdef G4VERBOSE
665 if(verboseLevel > 0)
666 { G4cout
667 << "#### G4VUserPhysicsList::BuildPhysicsTable() - BuildPhysicsTable("
668 << particle->GetParticleName() << ") skipped..." << G4endl; }
669#endif
670 return;
671 }
673 {
675 {
676 // fail to retrieve cut tables
677#ifdef G4VERBOSE
678 if(verboseLevel > 0)
679 {
680 G4cout << "G4VUserPhysicsList::BuildPhysicsTable "
681 << "Physics table can not be retrieved and will be calculated "
682 << G4endl;
683 }
684#endif
685 fRetrievePhysicsTable = false;
686 }
687 else
688 {
689#ifdef G4VERBOSE
690 if(verboseLevel > 2)
691 {
692 G4cout << "G4VUserPhysicsList::BuildPhysicsTable "
693 << " Retrieve Physics Table for " << particle->GetParticleName()
694 << G4endl;
695 }
696#endif
697 // Retrieve PhysicsTable from files for proccesses
699 }
700 }
701
702#ifdef G4VERBOSE
703 if(verboseLevel > 2)
704 {
705 G4cout << "G4VUserPhysicsList::BuildPhysicsTable "
706 << "Calculate Physics Table for " << particle->GetParticleName()
707 << G4endl;
708 }
709#endif
710 // Rebuild the physics tables for every process for this particle type
711 // if particle is not ShortLived
712 if(!particle->IsShortLived())
713 {
714 G4ProcessManager* pManager = particle->GetProcessManager();
715 if(pManager == nullptr)
716 {
717#ifdef G4VERBOSE
718 if(verboseLevel > 0)
719 {
720 G4cout << "G4VUserPhysicsList::BuildPhysicsTable "
721 << " : No Process Manager for " << particle->GetParticleName()
722 << G4endl;
723 G4cout << particle->GetParticleName()
724 << " should be created in your PhysicsList" << G4endl;
725 }
726#endif
727 G4Exception("G4VUserPhysicsList::BuildPhysicsTable", "Run0271",
728 FatalException, "No process manager");
729 return;
730 }
731
732 // Get processes from master thread;
733 G4ProcessManager* pManagerShadow = particle->GetMasterProcessManager();
734
735 G4ProcessVector* pVector = pManager->GetProcessList();
736 if(pVector == nullptr)
737 {
738#ifdef G4VERBOSE
739 if(verboseLevel > 0)
740 {
741 G4cout << "G4VUserPhysicsList::BuildPhysicsTable "
742 << " : No Process Vector for " << particle->GetParticleName()
743 << G4endl;
744 }
745#endif
746 G4Exception("G4VUserPhysicsList::BuildPhysicsTable", "Run0272",
747 FatalException, "No process Vector");
748 return;
749 }
750#ifdef G4VERBOSE
751 if(verboseLevel > 2)
752 {
753 G4cout << "G4VUserPhysicsList::BuildPhysicsTable %%%%%% "
754 << particle->GetParticleName() << G4endl;
755 G4cout << " ProcessManager : " << pManager
756 << " ProcessManagerShadow : " << pManagerShadow << G4endl;
757 for(std::size_t iv1 = 0; iv1 < pVector->size(); ++iv1)
758 {
759 G4cout << " " << iv1 << " - " << (*pVector)[iv1]->GetProcessName()
760 << G4endl;
761 }
762 G4cout << "--------------------------------------------------------------"
763 << G4endl;
764 G4ProcessVector* pVectorShadow = pManagerShadow->GetProcessList();
765
766 for(std::size_t iv2 = 0; iv2 < pVectorShadow->size(); ++iv2)
767 {
768 G4cout << " " << iv2 << " - "
769 << (*pVectorShadow)[iv2]->GetProcessName() << G4endl;
770 }
771 }
772#endif
773 for(std::size_t j = 0; j < pVector->size(); ++j)
774 {
775 // Andrea July 16th 2013 : migration to new interface...
776 // Infer if we are in a worker thread or master thread
777 // Master thread is the one in which the process manager
778 // and process manager shadow pointers are the same
779 if(pManagerShadow == pManager)
780 {
781 (*pVector)[j]->BuildPhysicsTable(*particle);
782 }
783 else
784 {
785 (*pVector)[j]->BuildWorkerPhysicsTable(*particle);
786 }
787
788 } // End loop on processes vector
789 } // End if short-lived
790}
@ FatalException
G4ProcessManager * GetProcessManager() const
G4VTrackingManager * GetTrackingManager() const
G4ProcessManager * GetMasterProcessManager() const
G4ProcessVector * GetProcessList() const
std::size_t size() const
virtual void RetrievePhysicsTable(G4ParticleDefinition *, const G4String &directory, G4bool ascii=false)

References directoryPhysicsTable, FatalException, fIsRestoredCutValues, fRetrievePhysicsTable, fStoredInAscii, G4cout, G4endl, G4Exception(), G4ParticleDefinition::GetMasterProcessManager(), G4ParticleDefinition::GetParticleName(), G4ProcessManager::GetProcessList(), G4ParticleDefinition::GetProcessManager(), G4ParticleDefinition::GetTrackingManager(), G4ParticleDefinition::IsShortLived(), RetrievePhysicsTable(), G4ProcessVector::size(), and verboseLevel.

◆ CheckParticleList()

void G4VUserPhysicsList::CheckParticleList ( )

Definition at line 1060 of file G4VUserPhysicsList.cc.

1061{
1063 {
1064 G4MT_thePLHelper->CheckParticleList();
1065 }
1066}

References fDisableCheckParticleList, and G4MT_thePLHelper.

Referenced by G4RunManagerKernel::InitializePhysics().

◆ Construct()

void G4VUserPhysicsList::Construct ( )
inline

Definition at line 319 of file G4VUserPhysicsList.hh.

320{
321 #ifdef G4VERBOSE
322 if(verboseLevel > 1)
323 G4cout << "G4VUserPhysicsList::Construct()" << G4endl;
324 #endif
325
327
329
330 #ifdef G4VERBOSE
331 if(verboseLevel > 1)
332 G4cout << "Construct processes " << G4endl;
333 #endif
335}
virtual void ConstructProcess()=0
G4bool IsMasterThread()
Definition: G4Threading.cc:124

References ConstructProcess(), G4cout, G4endl, G4PhysicsModelCatalog::Initialize(), InitializeProcessManager(), G4Threading::IsMasterThread(), and verboseLevel.

Referenced by G4RunManagerKernel::InitializePhysics().

◆ ConstructParticle()

virtual void G4VUserPhysicsList::ConstructParticle ( )
pure virtual

◆ ConstructProcess()

virtual void G4VUserPhysicsList::ConstructProcess ( )
pure virtual

◆ DisableCheckParticleList()

void G4VUserPhysicsList::DisableCheckParticleList ( )
inline

Definition at line 379 of file G4VUserPhysicsList.hh.

380{
382}

References fDisableCheckParticleList.

◆ DumpCutValuesTable()

void G4VUserPhysicsList::DumpCutValuesTable ( G4int  flag = 1)

◆ DumpCutValuesTableIfRequested()

void G4VUserPhysicsList::DumpCutValuesTableIfRequested ( )

◆ DumpList()

void G4VUserPhysicsList::DumpList ( ) const

Definition at line 895 of file G4VUserPhysicsList.cc.

896{
897 theParticleIterator->reset();
898 G4int idx = 0;
899 while((*theParticleIterator)())
900 {
901 G4ParticleDefinition* particle = theParticleIterator->value();
902 G4cout << particle->GetParticleName();
903 if((idx++ % 4) == 3)
904 {
905 G4cout << G4endl;
906 }
907 else
908 {
909 G4cout << ", ";
910 }
911 }
912 G4cout << G4endl;
913}
int G4int
Definition: G4Types.hh:85

References G4cout, G4endl, G4ParticleDefinition::GetParticleName(), and theParticleIterator.

Referenced by export_G4VUserPhysicsList(), and G4UserPhysicsListMessenger::SetNewValue().

◆ GetApplyCuts()

G4bool G4VUserPhysicsList::GetApplyCuts ( const G4String name) const

◆ GetCutValue()

G4double G4VUserPhysicsList::GetCutValue ( const G4String pname) const

Definition at line 421 of file G4VUserPhysicsList.cc.

422{
423 std::size_t nReg = (G4RegionStore::GetInstance())->size();
424 if(nReg == 0)
425 {
426#ifdef G4VERBOSE
427 if(verboseLevel > 0)
428 {
429 G4cout << "G4VUserPhysicsList::GetCutValue "
430 << " : No Default Region " << G4endl;
431 }
432#endif
433 G4Exception("G4VUserPhysicsList::GetCutValue", "Run0253", FatalException,
434 "No Default Region");
435 return -1. * mm;
436 }
437 G4Region* region =
438 G4RegionStore::GetInstance()->GetRegion("DefaultRegionForTheWorld", false);
439 return region->GetProductionCuts()->GetProductionCut(name);
440}
G4double GetProductionCut(G4int index) const
static G4RegionStore * GetInstance()
G4Region * GetRegion(const G4String &name, G4bool verbose=true) const
G4ProductionCuts * GetProductionCuts() const

References FatalException, G4cout, G4endl, G4Exception(), G4RegionStore::GetInstance(), G4ProductionCuts::GetProductionCut(), G4Region::GetProductionCuts(), G4RegionStore::GetRegion(), mm, G4InuclParticleNames::name(), and verboseLevel.

Referenced by SetCuts(), and G4UserPhysicsListMessenger::SetNewValue().

◆ GetDefaultCutValue()

G4double G4VUserPhysicsList::GetDefaultCutValue ( ) const
inline

Definition at line 337 of file G4VUserPhysicsList.hh.

338{
339 return defaultCutValue;
340}

References defaultCutValue.

Referenced by export_G4VUserPhysicsList(), and G4UserPhysicsListMessenger::GetCurrentValue().

◆ GetInstanceID()

G4int G4VUserPhysicsList::GetInstanceID ( ) const
inline

Definition at line 384 of file G4VUserPhysicsList.hh.

385{
386 return g4vuplInstanceID;
387}

References g4vuplInstanceID.

Referenced by G4VUserPhysicsList(), and operator=().

◆ GetParticleIterator()

G4ParticleTable::G4PTblDicIterator * G4VUserPhysicsList::GetParticleIterator ( ) const
protected

◆ GetPhysicsTableDirectory()

const G4String & G4VUserPhysicsList::GetPhysicsTableDirectory ( ) const
inline

◆ GetSubInstanceManager()

const G4VUPLManager & G4VUserPhysicsList::GetSubInstanceManager ( )
inlinestatic

◆ GetVerboseLevel()

G4int G4VUserPhysicsList::GetVerboseLevel ( ) const
inline

◆ InitializeProcessManager()

void G4VUserPhysicsList::InitializeProcessManager ( )
protected

Definition at line 217 of file G4VUserPhysicsList.cc.

218{
219 // Request lock for particle table access. Some changes are inside
220 // this critical region.
221#ifdef G4MULTITHREADED
222 G4MUTEXLOCK(&G4ParticleTable::particleTableMutex());
223 G4ParticleTable::lockCount()++;
224#endif
227
228 // loop over all particles in G4ParticleTable
229 theParticleIterator->reset();
230 while((*theParticleIterator)())
231 {
232 G4ParticleDefinition* particle = theParticleIterator->value();
233 G4ProcessManager* pmanager = particle->GetProcessManager();
234
235 if(pmanager == nullptr)
236 {
237 // create process manager if the particle does not have its own.
238 pmanager = new G4ProcessManager(particle);
239 particle->SetProcessManager(pmanager);
240 if(particle->GetMasterProcessManager() == nullptr)
241 particle->SetMasterProcessManager(pmanager);
242#ifdef G4VERBOSE
243 if(verboseLevel > 2)
244 {
245 G4cout << "G4VUserPhysicsList::InitializeProcessManager: creating "
246 "ProcessManager to "
247 << particle->GetParticleName() << G4endl;
248 }
249#endif
250 }
251 }
252
253 if(gion != nullptr)
254 {
255 G4ProcessManager* gionPM = gion->GetProcessManager();
256 // loop over all particles once again (this time, with all general ions)
257 theParticleIterator->reset(false);
258 while((*theParticleIterator)())
259 {
260 G4ParticleDefinition* particle = theParticleIterator->value();
261 if(particle->IsGeneralIon())
262 {
263 particle->SetProcessManager(gionPM);
264#ifdef G4VERBOSE
265 if(verboseLevel > 2)
266 {
267 G4cout << "G4VUserPhysicsList::InitializeProcessManager: copying "
268 "ProcessManager to "
269 << particle->GetParticleName() << G4endl;
270 }
271#endif
272 }
273 }
274 }
275
276 // release lock for particle table access
277#ifdef G4MULTITHREADED
278 G4MUTEXUNLOCK(&G4ParticleTable::particleTableMutex());
279#endif
280}
#define G4MUTEXLOCK(mutex)
Definition: G4Threading.hh:251
#define G4MUTEXUNLOCK(mutex)
Definition: G4Threading.hh:254
void SetMasterProcessManager(G4ProcessManager *aNewPM)
G4bool IsGeneralIon() const
void SetProcessManager(G4ProcessManager *aProcessManager)
G4ParticleDefinition * GetGenericIon() const

References G4cout, G4endl, G4MUTEXLOCK, G4MUTEXUNLOCK, G4ParticleTable::GetGenericIon(), G4ParticleDefinition::GetMasterProcessManager(), G4ParticleDefinition::GetParticleName(), G4ParticleTable::GetParticleTable(), G4ParticleDefinition::GetProcessManager(), G4ParticleDefinition::IsGeneralIon(), G4ParticleDefinition::SetMasterProcessManager(), G4ParticleDefinition::SetProcessManager(), theParticleIterator, and verboseLevel.

Referenced by Construct().

◆ InitializeWorker()

void G4VUserPhysicsList::InitializeWorker ( )
virtual

Definition at line 117 of file G4VUserPhysicsList.cc.

118{
119 // Remember messengers are per-thread, so this needs to be done by each
120 // worker and due to the presence of "this" cannot be done in
121 // G4VUPLData::initialize()
123}

References G4MT_theMessenger.

Referenced by G4WorkerRunManager::SetUserInitialization().

◆ IsPhysicsTableRetrieved()

G4bool G4VUserPhysicsList::IsPhysicsTableRetrieved ( ) const
inline

◆ IsStoredInAscii()

G4bool G4VUserPhysicsList::IsStoredInAscii ( ) const
inline

Definition at line 352 of file G4VUserPhysicsList.hh.

353{
354 return fStoredInAscii;
355}

References fStoredInAscii.

Referenced by export_G4VUserPhysicsList(), and G4UserPhysicsListMessenger::GetCurrentValue().

◆ operator=()

G4VUserPhysicsList & G4VUserPhysicsList::operator= ( const G4VUserPhysicsList right)

◆ PreparePhysicsTable()

void G4VUserPhysicsList::PreparePhysicsTable ( G4ParticleDefinition particle)

Definition at line 793 of file G4VUserPhysicsList.cc.

794{
795 if (auto *trackingManager = particle->GetTrackingManager())
796 {
797 trackingManager->PreparePhysicsTable(*particle);
798 return;
799 }
800
801 if(!(particle->GetMasterProcessManager()))
802 {
803 return;
804 }
805 // Prepare the physics tables for every process for this particle type
806 // if particle is not ShortLived
807 if(!particle->IsShortLived())
808 {
809 G4ProcessManager* pManager = particle->GetProcessManager();
810 if(pManager == nullptr)
811 {
812#ifdef G4VERBOSE
813 if(verboseLevel > 0)
814 {
815 G4cout << "G4VUserPhysicsList::PreparePhysicsTable "
816 << ": No Process Manager for " << particle->GetParticleName()
817 << G4endl;
818 G4cout << particle->GetParticleName()
819 << " should be created in your PhysicsList" << G4endl;
820 }
821#endif
822 G4Exception("G4VUserPhysicsList::PreparePhysicsTable", "Run0273",
823 FatalException, "No process manager");
824 return;
825 }
826
827 // Get processes from master thread
828 G4ProcessManager* pManagerShadow = particle->GetMasterProcessManager();
829
830 G4ProcessVector* pVector = pManager->GetProcessList();
831 if(pVector == nullptr)
832 {
833#ifdef G4VERBOSE
834 if(verboseLevel > 0)
835 {
836 G4cout << "G4VUserPhysicsList::PreparePhysicsTable "
837 << ": No Process Vector for " << particle->GetParticleName()
838 << G4endl;
839 }
840#endif
841 G4Exception("G4VUserPhysicsList::PreparePhysicsTable", "Run0274",
842 FatalException, "No process Vector");
843 return;
844 }
845 for(std::size_t j = 0; j < pVector->size(); ++j)
846 {
847 // Andrea July 16th 2013 : migration to new interface...
848 // Infer if we are in a worker thread or master thread
849 // Master thread is the one in which the process manager
850 // and process manager shadow pointers are the same
851 if(pManagerShadow == pManager)
852 {
853 (*pVector)[j]->PreparePhysicsTable(*particle);
854 }
855 else
856 {
857 (*pVector)[j]->PrepareWorkerPhysicsTable(*particle);
858 }
859 } // End loop on processes vector
860 } // End if pn ShortLived
861}

References FatalException, G4cout, G4endl, G4Exception(), G4ParticleDefinition::GetMasterProcessManager(), G4ParticleDefinition::GetParticleName(), G4ProcessManager::GetProcessList(), G4ParticleDefinition::GetProcessManager(), G4ParticleDefinition::GetTrackingManager(), G4ParticleDefinition::IsShortLived(), G4ProcessVector::size(), and verboseLevel.

Referenced by BuildPhysicsTable(), and G4UserPhysicsListMessenger::SetNewValue().

◆ RegisterProcess()

G4bool G4VUserPhysicsList::RegisterProcess ( G4VProcess process,
G4ParticleDefinition particle 
)
protected

Definition at line 1081 of file G4VUserPhysicsList.cc.

1083{
1084 return G4MT_thePLHelper->RegisterProcess(process, particle);
1085}

References G4MT_thePLHelper.

◆ RemoveProcessManager()

void G4VUserPhysicsList::RemoveProcessManager ( )

Definition at line 283 of file G4VUserPhysicsList.cc.

284{
285 // Request lock for particle table access. Some changes are inside
286 // this critical region.
287#ifdef G4MULTITHREADED
288 G4MUTEXLOCK(&G4ParticleTable::particleTableMutex());
289 G4ParticleTable::lockCount()++;
290#endif
291
292 // loop over all particles in G4ParticleTable
293 theParticleIterator->reset();
294 while((*theParticleIterator)())
295 {
296 G4ParticleDefinition* particle = theParticleIterator->value();
297 if(particle->GetInstanceID() <
299 {
300 if(particle->GetParticleSubType() != "generic" ||
301 particle->GetParticleName() == "GenericIon")
302 {
303 G4ProcessManager* pmanager = particle->GetProcessManager();
304 if(pmanager != nullptr)
305 delete pmanager;
306#ifdef G4VERBOSE
307 if(verboseLevel > 2)
308 {
309 G4cout << "G4VUserPhysicsList::RemoveProcessManager: ";
310 G4cout << "remove ProcessManager from ";
311 G4cout << particle->GetParticleName() << G4endl;
312 }
313#endif
314 }
315 particle->SetProcessManager(nullptr);
316 }
317 }
318
319 // release lock for particle table access
320#ifdef G4MULTITHREADED
321 G4MUTEXUNLOCK(&G4ParticleTable::particleTableMutex());
322#endif
323}
static G4PART_DLL G4int & slavetotalspace()
G4int GetInstanceID() const
const G4String & GetParticleSubType() const

References G4cout, G4endl, G4MUTEXLOCK, G4MUTEXUNLOCK, G4ParticleDefinition::GetInstanceID(), G4ParticleDefinition::GetParticleName(), G4ParticleDefinition::GetParticleSubType(), G4ParticleDefinition::GetProcessManager(), G4ParticleDefinition::SetProcessManager(), G4PDefManager::slavetotalspace(), theParticleIterator, and verboseLevel.

Referenced by TerminateWorker(), and ~G4VUserPhysicsList().

◆ RemoveTrackingManager()

void G4VUserPhysicsList::RemoveTrackingManager ( )

Definition at line 326 of file G4VUserPhysicsList.cc.

327{
328 // One tracking manager may be registered for multiple particles, make sure
329 // to delete every object only once.
330 std::unordered_set<G4VTrackingManager *> trackingManagers;
331
332 // loop over all particles in G4ParticleTable
333 theParticleIterator->reset();
334 while((*theParticleIterator)())
335 {
336 G4ParticleDefinition* particle = theParticleIterator->value();
337 if (auto *trackingManager = particle->GetTrackingManager())
338 {
339#ifdef G4VERBOSE
340 if(verboseLevel > 2)
341 {
342 G4cout << "G4VUserPhysicsList::RemoveTrackingManager: ";
343 G4cout << "remove TrackingManager from ";
344 G4cout << particle->GetParticleName() << G4endl;
345 }
346#endif
347 trackingManagers.insert(trackingManager);
348 particle->SetTrackingManager(nullptr);
349 }
350 }
351
352 for (G4VTrackingManager *tm : trackingManagers)
353 {
354 delete tm;
355 }
356}
void SetTrackingManager(G4VTrackingManager *aTrackingManager)

References G4cout, G4endl, G4ParticleDefinition::GetParticleName(), G4ParticleDefinition::GetTrackingManager(), G4ParticleDefinition::SetTrackingManager(), theParticleIterator, G4InuclParticleNames::tm, and verboseLevel.

Referenced by TerminateWorker(), and ~G4VUserPhysicsList().

◆ ResetPhysicsTableRetrieved()

void G4VUserPhysicsList::ResetPhysicsTableRetrieved ( )
inline

◆ ResetStoredInAscii()

void G4VUserPhysicsList::ResetStoredInAscii ( )
inline

Definition at line 374 of file G4VUserPhysicsList.hh.

375{
376 fStoredInAscii = false;
377}

References fStoredInAscii.

Referenced by export_G4VUserPhysicsList(), and G4UserPhysicsListMessenger::SetNewValue().

◆ RetrievePhysicsTable()

void G4VUserPhysicsList::RetrievePhysicsTable ( G4ParticleDefinition particle,
const G4String directory,
G4bool  ascii = false 
)
protectedvirtual

Definition at line 996 of file G4VUserPhysicsList.cc.

999{
1000 G4bool success[100];
1001 // Retrieve physics tables for every process for this particle type
1002 G4ProcessVector* pVector = (particle->GetProcessManager())->GetProcessList();
1003 for(std::size_t j = 0; j < pVector->size(); ++j)
1004 {
1005 success[j] =
1006 (*pVector)[j]->RetrievePhysicsTable(particle, directory, ascii);
1007
1008 if(!success[j])
1009 {
1010#ifdef G4VERBOSE
1011 if(verboseLevel > 2)
1012 {
1013 G4cout << "G4VUserPhysicsList::RetrievePhysicsTable "
1014 << " Fail to retrieve Physics Table for "
1015 << (*pVector)[j]->GetProcessName() << G4endl;
1016 G4cout << "Calculate Physics Table for " << particle->GetParticleName()
1017 << G4endl;
1018 }
1019#endif
1020 (*pVector)[j]->BuildPhysicsTable(*particle);
1021 }
1022 }
1023 for(std::size_t j = 0; j < pVector->size(); ++j)
1024 {
1025 // temporary addition to make the integral schema
1026 if(!success[j])
1027 BuildIntegralPhysicsTable((*pVector)[j], particle);
1028 }
1029}
bool G4bool
Definition: G4Types.hh:86
void BuildIntegralPhysicsTable(G4VProcess *, G4ParticleDefinition *)

References BuildIntegralPhysicsTable(), G4cout, G4endl, G4ParticleDefinition::GetParticleName(), G4ParticleDefinition::GetProcessManager(), G4ProcessVector::size(), and verboseLevel.

Referenced by BuildPhysicsTable().

◆ SetApplyCuts()

void G4VUserPhysicsList::SetApplyCuts ( G4bool  value,
const G4String name 
)

Definition at line 1032 of file G4VUserPhysicsList.cc.

1033{
1034#ifdef G4VERBOSE
1035 if(verboseLevel > 2)
1036 {
1037 G4cout << "G4VUserPhysicsList::SetApplyCuts for " << name << G4endl;
1038 }
1039#endif
1040 if(name == "all")
1041 {
1046 }
1047 else
1048 {
1050 }
1051}

References G4ParticleTable::FindParticle(), G4cout, G4endl, G4InuclParticleNames::name(), G4ParticleDefinition::SetApplyCutsFlag(), theParticleTable, and verboseLevel.

Referenced by G4UserPhysicsListMessenger::SetNewValue().

◆ SetCuts()

void G4VUserPhysicsList::SetCuts ( )
virtual

Reimplemented in PhysicsList, pyG4VModularPhysicsList::CB_G4VModularPhysicsList, pyG4VUserPhysicsList::CB_G4VUserPhysicsList, PhysicsListEMstd, ExN01PhysicsList, QPhysicsList, G4ErrorPhysicsList, FTFP_BERT_HP, LBE, and QGSP_BERT_HP.

Definition at line 359 of file G4VUserPhysicsList.cc.

360{
362 {
364 }
365
366#ifdef G4VERBOSE
367 if(verboseLevel > 1)
368 {
369 G4cout << "G4VUserPhysicsList::SetCuts: " << G4endl;
370 G4cout << "Cut for gamma: " << GetCutValue("gamma") / mm << "[mm]"
371 << G4endl;
372 G4cout << "Cut for e-: " << GetCutValue("e-") / mm << "[mm]" << G4endl;
373 G4cout << "Cut for e+: " << GetCutValue("e+") / mm << "[mm]" << G4endl;
374 G4cout << "Cut for proton: " << GetCutValue("proton") / mm << "[mm]"
375 << G4endl;
376 }
377#endif
378
379 // dump Cut values if verboseLevel==3
380 if(verboseLevel > 2)
381 {
383 }
384}
G4double GetCutValue(const G4String &pname) const
void SetDefaultCutValue(G4double newCutValue)
void DumpCutValuesTable(G4int flag=1)

References defaultCutValue, DumpCutValuesTable(), G4cout, G4endl, GetCutValue(), isSetDefaultCutValue, mm, SetDefaultCutValue(), and verboseLevel.

Referenced by export_G4VModularPhysicsList(), export_G4VUserPhysicsList(), G4RunManagerKernel::InitializePhysics(), SetCutsWithDefault(), and G4UserPhysicsListMessenger::SetNewValue().

◆ SetCutsForRegion()

void G4VUserPhysicsList::SetCutsForRegion ( G4double  aCut,
const G4String rname 
)

Definition at line 478 of file G4VUserPhysicsList.cc.

479{
480 // set cut values for gamma at first and for e- and e+
481 SetCutValue(aCut, "gamma", rname);
482 SetCutValue(aCut, "e-", rname);
483 SetCutValue(aCut, "e+", rname);
484 SetCutValue(aCut, "proton", rname);
485}
void SetCutValue(G4double aCut, const G4String &pname)

References SetCutValue().

Referenced by export_G4VUserPhysicsList(), and G4UserPhysicsListMessenger::SetNewValue().

◆ SetCutsWithDefault()

void G4VUserPhysicsList::SetCutsWithDefault ( )

◆ SetCutValue() [1/2]

void G4VUserPhysicsList::SetCutValue ( G4double  aCut,
const G4String pname 
)

◆ SetCutValue() [2/2]

void G4VUserPhysicsList::SetCutValue ( G4double  aCut,
const G4String pname,
const G4String rname 
)

Definition at line 449 of file G4VUserPhysicsList.cc.

451{
453 if(region != nullptr)
454 {
455 // set cut value
456 SetParticleCuts(aCut, pname, region);
457 }
458 else
459 {
460#ifdef G4VERBOSE
461 if(verboseLevel > 0)
462 {
463 G4cout << "G4VUserPhysicsList::SetCutValue "
464 << " : No Region of " << rname << G4endl;
465 }
466#endif
467 }
468}
string pname
Definition: eplot.py:33

References G4cout, G4endl, G4RegionStore::GetInstance(), G4RegionStore::GetRegion(), eplot::pname, SetParticleCuts(), and verboseLevel.

◆ SetDefaultCutValue()

void G4VUserPhysicsList::SetDefaultCutValue ( G4double  newCutValue)

Definition at line 387 of file G4VUserPhysicsList.cc.

388{
389 if(value < 0.0)
390 {
391#ifdef G4VERBOSE
392 if(verboseLevel > 0)
393 {
394 G4cout << "G4VUserPhysicsList::SetDefaultCutValue: negative cut values"
395 << " :" << value / mm << "[mm]" << G4endl;
396 }
397#endif
398 return;
399 }
400
401 defaultCutValue = value;
403
404 // set cut values for gamma at first and for e- and e+
408 SetCutValue(defaultCutValue, "proton");
409
410#ifdef G4VERBOSE
411 if(verboseLevel > 1)
412 {
413 G4cout << "G4VUserPhysicsList::SetDefaultCutValue:"
414 << "default cut value is changed to :" << defaultCutValue / mm
415 << "[mm]" << G4endl;
416 }
417#endif
418}

References defaultCutValue, G4cout, G4endl, isSetDefaultCutValue, mm, SetCutValue(), and verboseLevel.

Referenced by export_G4VUserPhysicsList(), SetCuts(), SetCutsWithDefault(), G4UserPhysicsListMessenger::SetNewValue(), and SetParticleCuts().

◆ SetParticleCuts() [1/2]

void G4VUserPhysicsList::SetParticleCuts ( G4double  cut,
const G4String particleName,
G4Region region = nullptr 
)

Definition at line 496 of file G4VUserPhysicsList.cc.

499{
500 if(cut < 0.0)
501 {
502#ifdef G4VERBOSE
503 if(verboseLevel > 0)
504 {
505 G4cout << "G4VUserPhysicsList::SetParticleCuts: negative cut values"
506 << " :" << cut / mm << "[mm]"
507 << " for " << particleName << G4endl;
508 }
509#endif
510 return;
511 }
512
513 G4Region* world_region =
514 G4RegionStore::GetInstance()->GetRegion("DefaultRegionForTheWorld", false);
515 if(region == nullptr)
516 {
517 std::size_t nReg = G4RegionStore::GetInstance()->size();
518 if(nReg == 0)
519 {
520#ifdef G4VERBOSE
521 if(verboseLevel > 0)
522 {
523 G4cout << "G4VUserPhysicsList::SetParticleCuts "
524 << " : No Default Region " << G4endl;
525 }
526#endif
527 G4Exception("G4VUserPhysicsList::SetParticleCuts ", "Run0254",
528 FatalException, "No Default Region");
529 return;
530 }
531 region = world_region;
532 }
533
535 {
537 }
538
539 G4ProductionCuts* pcuts = region->GetProductionCuts();
540 if(region != world_region &&
542 ->GetDefaultProductionCuts())
543 { // This region had no unique cuts yet but shares the default cuts.
544 // Need to create a new object before setting the value.
545 pcuts =
547 ->GetDefaultProductionCuts()));
548 region->SetProductionCuts(pcuts);
549 }
550 pcuts->SetProductionCut(cut, particleName);
551#ifdef G4VERBOSE
552 if(verboseLevel > 2)
553 {
554 G4cout << "G4VUserPhysicsList::SetParticleCuts: "
555 << " :" << cut / mm << "[mm]"
556 << " for " << particleName << G4endl;
557 }
558#endif
559}
void SetProductionCut(G4double cut, G4int index=-1)
void SetProductionCuts(G4ProductionCuts *cut)

References defaultCutValue, FatalException, G4cout, G4endl, G4Exception(), G4RegionStore::GetInstance(), G4Region::GetProductionCuts(), G4ProductionCutsTable::GetProductionCutsTable(), G4RegionStore::GetRegion(), isSetDefaultCutValue, mm, SetDefaultCutValue(), G4ProductionCuts::SetProductionCut(), G4Region::SetProductionCuts(), and verboseLevel.

◆ SetParticleCuts() [2/2]

void G4VUserPhysicsList::SetParticleCuts ( G4double  cut,
G4ParticleDefinition particle,
G4Region region = nullptr 
)

◆ SetPhysicsTableRetrieved()

void G4VUserPhysicsList::SetPhysicsTableRetrieved ( const G4String directory = "")

Definition at line 984 of file G4VUserPhysicsList.cc.

985{
987 if(!directory.empty())
988 {
989 directoryPhysicsTable = directory;
990 }
992 fIsRestoredCutValues = false;
993}

References directoryPhysicsTable, fIsCheckedForRetrievePhysicsTable, fIsRestoredCutValues, and fRetrievePhysicsTable.

Referenced by G4UserPhysicsListMessenger::SetNewValue().

◆ SetStoredInAscii()

void G4VUserPhysicsList::SetStoredInAscii ( )
inline

Definition at line 362 of file G4VUserPhysicsList.hh.

363{
364 fStoredInAscii = true;
365}

References fStoredInAscii.

Referenced by export_G4VUserPhysicsList(), and G4UserPhysicsListMessenger::SetNewValue().

◆ SetVerboseLevel()

void G4VUserPhysicsList::SetVerboseLevel ( G4int  value)

Definition at line 1095 of file G4VUserPhysicsList.cc.

1096{
1097 verboseLevel = value;
1098 // set verboseLevel for G4ProductionCutsTable same as one for
1099 // G4VUserPhysicsList:
1101
1102 G4MT_thePLHelper->SetVerboseLevel(verboseLevel);
1103
1104#ifdef G4VERBOSE
1105 if(verboseLevel > 1)
1106 {
1107 G4cout << "G4VUserPhysicsList::SetVerboseLevel :"
1108 << " Verbose level is set to " << verboseLevel << G4endl;
1109 }
1110#endif
1111}
void SetVerboseLevel(G4int value)

References fCutsTable, G4cout, G4endl, G4MT_thePLHelper, G4ProductionCutsTable::SetVerboseLevel(), and verboseLevel.

Referenced by export_G4VUserPhysicsList(), ExN01PhysicsList::SetCuts(), and G4UserPhysicsListMessenger::SetNewValue().

◆ StorePhysicsTable()

G4bool G4VUserPhysicsList::StorePhysicsTable ( const G4String directory = ".")

Definition at line 931 of file G4VUserPhysicsList.cc.

932{
933 G4bool ascii = fStoredInAscii;
934 G4String dir = directory;
935 if(dir.empty())
937 else
939
940 // store CutsTable info
941 if(!fCutsTable->StoreCutsTable(dir, ascii))
942 {
943 G4Exception("G4VUserPhysicsList::StorePhysicsTable", "Run0281", JustWarning,
944 "Fail to store Cut Table");
945 return false;
946 }
947#ifdef G4VERBOSE
948 if(verboseLevel > 2)
949 {
950 G4cout << "G4VUserPhysicsList::StorePhysicsTable "
951 << " Store material and cut values successfully" << G4endl;
952 }
953#endif
954
955 G4bool success = true;
956
957 // loop over all particles in G4ParticleTable
958 theParticleIterator->reset();
959 while((*theParticleIterator)())
960 {
961 G4ParticleDefinition* particle = theParticleIterator->value();
962 // Store physics tables for every process for this particle type
963 G4ProcessVector* pVector =
964 (particle->GetProcessManager())->GetProcessList();
965 for(std::size_t j = 0; j < pVector->size(); ++j)
966 {
967 if(!(*pVector)[j]->StorePhysicsTable(particle, dir, ascii))
968 {
969 G4String comment = "Fail to store physics table for ";
970 comment += (*pVector)[j]->GetProcessName();
971 comment += "(" + particle->GetParticleName() + ")";
972 G4Exception("G4VUserPhysicsList::StorePhysicsTable", "Run0282",
973 JustWarning, comment);
974 success = false;
975 }
976 }
977 // end loop over processes
978 }
979 // end loop over particles
980 return success;
981}
G4bool StoreCutsTable(const G4String &directory, G4bool ascii=false)

References directoryPhysicsTable, fCutsTable, fStoredInAscii, G4cout, G4endl, G4Exception(), G4ParticleDefinition::GetParticleName(), G4ParticleDefinition::GetProcessManager(), JustWarning, G4ProcessVector::size(), G4ProductionCutsTable::StoreCutsTable(), theParticleIterator, and verboseLevel.

Referenced by export_G4VUserPhysicsList(), and G4UserPhysicsListMessenger::SetNewValue().

◆ TerminateWorker()

void G4VUserPhysicsList::TerminateWorker ( )
virtual

◆ UseCoupledTransportation()

void G4VUserPhysicsList::UseCoupledTransportation ( G4bool  vl = true)

Definition at line 1075 of file G4VUserPhysicsList.cc.

1076{
1077 G4MT_thePLHelper->UseCoupledTransportation(vl);
1078}

References G4MT_thePLHelper.

Referenced by G4RunManagerKernel::InitializePhysics().

Field Documentation

◆ defaultCutValue

G4double G4VUserPhysicsList::defaultCutValue = 1.0
protected

◆ directoryPhysicsTable

G4String G4VUserPhysicsList::directoryPhysicsTable = "."
protected

◆ fCutsTable

G4ProductionCutsTable* G4VUserPhysicsList::fCutsTable = nullptr
protected

◆ fDisableCheckParticleList

G4bool G4VUserPhysicsList::fDisableCheckParticleList = false
protected

◆ fIsCheckedForRetrievePhysicsTable

G4bool G4VUserPhysicsList::fIsCheckedForRetrievePhysicsTable = false
protected

◆ fIsRestoredCutValues

G4bool G4VUserPhysicsList::fIsRestoredCutValues = false
protected

◆ fRetrievePhysicsTable

G4bool G4VUserPhysicsList::fRetrievePhysicsTable = false
protected

◆ fStoredInAscii

G4bool G4VUserPhysicsList::fStoredInAscii = true
protected

◆ g4vuplInstanceID

G4int G4VUserPhysicsList::g4vuplInstanceID = 0
protected

◆ isSetDefaultCutValue

G4bool G4VUserPhysicsList::isSetDefaultCutValue = false
protected

◆ subInstanceManager

G4VUPLManager G4VUserPhysicsList::subInstanceManager
staticprotected

◆ theParticleTable

G4ParticleTable* G4VUserPhysicsList::theParticleTable = nullptr
protected

◆ verboseLevel

G4int G4VUserPhysicsList::verboseLevel = 1
protected

The documentation for this class was generated from the following files: